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Abstract 

We develop a semi-quantitative theory of electron pairing and resulting superconduc- 
tivity in bulk "poor conductors" in which Fermi energy Ep is located in the region of 
localized states not so far from the Anderson mobility edge Ec- We assume attractive 
interaction between electrons near the Fermi surface. We review the existing theories 
and experimental data and argue that a large class of disordered films is described by 
this model. 

Our theoretical analysis is based on analytical treatment of pairing correlations, de- 
scribed in the basis of the exact single-particle eigenstates of the 3D Anderson model, 
which we combine with numerical data on eigenfunction correlations. Fractal nature of 
critical wavefunction's correlations is shown to be crucial for the physics of these systems. 

We identify three distinct phases: 'critical' superconductive state formed at Ep = E^ 
superconducting state with a strong pseudogap, realized due to pairing of weakly localized 
electrons and insulating state realized at Ep still deeper inside localized band. The 
'critical' superconducting phase is characterized by the enhancement of the transition 
temperature with respect to BCS result, by the inhomogeneous spatial distribution of 
superconductive order parameter and local density of states. The major new feature of 
the pseudo-gaped state is the presence of two independent energy scales: superconducting 
gap A, that is due to many-body correlations and a new "pseudogap" energy scale Ap 
which characterizes typical binding energy of localized electron pairs and leads to the 
insulating behavior of the resistivity as a function of temperature above superconductive 
Tc- Two gap nature of the pseudogapped superconductor is shown to lead to specific 
features seen in scanning tunneling spectroscopy and point-contact Andreev spectroscopy. 
We predict that pseudogaped superconducting state demonstrates anomalous behavior of 
the optical spectral weight. The insulating state is realized due to presence of local pairing 
gap but without superconducting correlations; it is characterized by a hard insulating 
gap in the density of single electrons and by purely activated low-temperature resistivity 
In i?(T) ~ l/T. 

Based on these results we propose a new "pseudospin" scenario of superconductor- 
insulator transition and argue that it is realized in a particular class of disordered super- 
conducting films. We conclude by the discussion of the experimental predictions of the 
theory and the theoretical issues that remain unsolved. 
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1. Introduction 

The purpose of this paper is to develop the semi-quantitative extension of the BCS 
theory of superconductivity that describes strongly disordered conductors which normal 
state is a weak Anderson insulator [l[ or a very poor metal. The paper focuses on the case 
of "uniformly disordered" materials, which do not contain morphological structures such 
as grains coupled by tunnel junctions. Below in this section we briefly review the existing 
theoretical models of the superconductor-insulator transition (SIT), compare them with 
the results of the experimental studies of the uniformly disordered films and choose 
the appropriate model for the superconductor-insulator transition in these materials. 
The conclusion of this introductory part is that this quantum transition in uniformly 
disordered films can be described by BCS pairing of electrons which single-particle states 
are close to the mobility edge of Anderson localization [l[. Because BCS pairing is most 
relevant for electrons close to the Fermi surface, the transition occurs when Fermi-level 
Ep is located in the region of localized single-electron states but close to the mobility 
edge. In the vicinity of the transition localization length is longer than typical distance 
between carriers while the relevant single electron states have the statistical properties 
of " critical wavefunctions" . Section [5] formulates the model in more detail and discusses 
the issue of wavefunction fractality, which turns out to be very important for the theory 
of superconductor-insulator transition developed in this work. 



3 



The important difference between our approach and many other works on superconductor- 
insulator transition is neglect of the effects of Coulomb interaction but consistent treat- 
ment of moderately strong disorder. In this respect our work is the extension of the 
approaches developed originally by Ma and Lee Kapitulnik and Kotliar 01, Bu- 
laevskii and Sadovskii [ij, and more recently by Ghosal, Trivcdi and Randcria |5| that 
have considered competition of superconducting pairing and Anderson localization with- 
out explicit account for Coulomb interaction. We give detailed arguments which justify 
applicability of our approach to disordered films such as amorphous InOa; and TiN on 
both phenomenological and microscopic levels below in subsection ll.2l Briefly, one should 
distinguish two possible effects of the Coulomb interaction: suppression of paring inter- 
action between individual fermions which occurs at short scales and enhancement of the 
phase fluctuations of the order parameter at large scales. The first would lead to a gap 
suppression in a direct contradiction with data while the second would lead to the phe- 
nomenology similar to that of Josephson junction arrays which display markedly different 
behavior. Our theory can be also applied to the cold atoms in optical and magnetic lat- 
tices with a controlled disorder. In such systems interaction is always attractive and 
effectively short-range. Moreover, it is tunable by magnetic field due to the Feshbach 
resonance Q which would allow to test directly the detailed predictions of the developed 
theory. 

The important consequence of the wave function fractality is the formation of the 
strongly bound electron pairs which survive deep in the insulating regime. In this 
situation the pairing interaction reduces the mobility of individual electons leading to 
" superconductivity- induced" insulator. We discuss this behavior in section [3] Theory of 
Cooper pairing of electrons populating critical fractal states is developed in section 2] 
Here we develop three different approximations for the computation of the superconduct- 
ing transition temperature and other properties. Section [5] gives the physical properties 
of the superconductor in this regime. We show that three different approximations devel- 
oped in section |4] agree with each other and predict parametrically strong enhancement 
of Tc with respect to its value given by the " Anderson theorem" . Another distinguishing 
feature of the emerging superconducting state is extremely strong spatial inhomogeneity 
of superconducting order parameter and the presence of a well-defined global T^. 

Section [5] presents the theory of superconductivity coexisting with a strong pseudo- 
gap. Here we show that superconductivity survives when Ep is located much deeper in 
the localized band than it was previously expected. In this regime the superconductiv- 
ity develops against the background of the pseudogap and, thus, is characterized by a 
number of unconventional properties. We present the specific results for the tunneling 
density of states and tunneling conductance, Andreev point-contact conductance, and 
spectral weight of high-frequency conductivity in this pseudogap superconductive state. 
Finally, section [7] reviews the main results and discusses a number of open problems. 
The Appendices present technical details of virial expansion that was used as the one of 
methods for the determination of Tc- 

1.1. Theoretical models 

Quantum phase transitions from superconducting to insulating state in disordered 
conductors and artificial structures were studied intensively since mid-1980s, for review 
see e.g. 0, H 0- A number of theoretical models describing such transitions were 
proposed and studied but a fully coherent theoretical picture of this phenomena has 
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not been established. We believe that such transitions might be driven by different 
mechanisms in different materials and thus belong to different universality classes which 
should be described by different models. Below in this section wc briefly review the 
alternative mechanisms and the theoretical methods employed for their description (see 
also short review [l^l)- To avoid confusion, we note that the first and the second of 
the mechanisms described below will not be studied in the main part of our paper. We 
describe them in some detail mostly because we need these details in order to argue below 
in this section that they are not relevant for the superconductor-insulator transition in 
homogeneous amorphous films. 

1.1.1. Coulomb blockade versus Superconductivity in J osephson junction arrays. 

Macroscopic conductor that looks homogeneous at a macroscopic scale might be in 
fact composed of small grains or islands of a good superconducting metal with tran- 
sition temperature Tcq. These grains are coupled to each other via low-transparency 
insulating tunnel barriers 0, 0, characterized by dimensionless tunnel conductances 
Gij = h/{2e)^Rij. The transition is due to the competition of the charging and Joseph- 
son energies. Exactly the same physics is realized in the artificial Josephson junction 
arrays, the only difference between the inhomogeneous films and artificial array is that 
the number of neighbors and of plaquette areas in the former are random. The simplest 
model describing this physics is given by effective Hamiltonian written in terms of phases 

and charges Qj assigned to grains: 



where Nj = Qj/2e is the number of Cooper pairs on the j-th grain, Cij is the matrix 
of mutual capacitances, and Ey is the Josephson coupling energy, correspondingly. The 
model is often further simplified by assuming that Ey = Ej is a non-zero constant only 
for nearest-neighboring grains while the capacitance matrix is diagonal, Cij = CoSij . This 
model neglects the effects of the quasiparticles which might be (sometimes) justified for 
small grains at low temperatures T <§; Tco, due to exponentially small density of normal 
electrons. The key parameter of the problem is the energy ratio x = Ej/Ec where 
Ec ~ (2e)^/2Co is the Coulomb charging energy due to 2e charge transfer. As was 
shown in the paper [ll| the ground-state of the Hamiltonian ([T]) is insulating at x ^ 1 
and superconducting at a; 3> 1, thus a phase transition(s) takes place at a; ^ 1. An 
essence of this phase transition is the Mott-Hubbard localization of Cooper pairs, taking 
place when tunneling matrix element of a pair {Ej) is much less than on-site repulsion 



However, the model (HJ is unlikely to describe correctly the physics of superconductor- 
insulator transition in Josephson arrays or granular materials, especially in its simplified 
version with diagonal capacitance matrix. It has two important deficiencies. 

First, realistic Josephson junction arrays and disordered films are poorly described 
by the model of diagonal capacitance matrix Cij = CgSij because normally the charging 
effects are controlled by capacitances of junctions C ^ Cp, not by the ground capaci- 
tances of the islands (see Ref. Q ) . It is in fact impossible to have a capacitance matrix 
dominated by the ground capacitance in the arrays which dimensionless normal state 
conductance is G = /i/[(2e)^i?T] ^ 1 because in these arrays the capacitance of the 




(1) 



Ec- 
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junctions cannot be small. The reason for this is that apart from purely geometrical 
contribution Cgcom = ^T^S/d, junction capacitance C = -\- C™'^ contains additional 

induced term Cind = ^Ce^/A, (this expression is valid at T ^ A) This induced con- 



tribution is due to virtual electron transition across the gap [l2|, [l3|. As a result, the 
charging energy can not be made arbitrary large (equivalently capacitance cannot be 
small): Ec < 32A/3G'. Josephson energy of the symmetric junction at low temperatures 
is Ej = GA/2, thus the condition Eq > Ej cannot be realized at large G. Furthermore, 



at temperatures above the parity effect threshold T* , (see jl4l | ) an additional contribution 
to the screening of Coulomb interaction between Cooper pairs comes from single-electron 
tunneling. Thus, in all cases the effect of capacitance renormalization is that the ratio 
X = Ej /Ec in granular arrays is controlled by the dimensionless conductance G in such 
a way that Coulomb effects are always weak at G > 1. Ground capacitance larger than 
the junction capacitance thus implies that the transition into the insulating state would 
occur in the arrays characterized by very small G ^ 1. Such behavior was never observed 
experimentally in Josephson arrays (see section rO)) . 

More realistic model involves the capacitance matrix that is dominated by the junc- 
tions capacitances with a small contribution from the ground capacitance of each grain. 
In this case the arguments of preceding paragraph show that the transition between the 
insulating and superconducting state should occur at Gc ~ l.fljj Moreover, deep in 
the insulating phase the electrostatic interaction between the charges in 2D Josephson 
array becomes logarithmic in distance, similar to the one between the vortices in the 
superconducting phase. Assumption of the full duality between vortices in the super- 
conducting state and charges in the insulating state allows one to make a number of 
predictions. [l6, 17 1 For instance, because the current of vortices generates voltage while 



the current of pairs implies the electrical current one expects that sup erconducting 



insulator transition is characterized by the universal value of Gc = 1. 17 

Unfortunately, despite a significant experimental effort the universal value of the 
resistance was never experimentally confirmed for Josephson arrays and for most dis- 
ordered films. Moreover, in non-zero magnetic field the Josephson arrays often show a 
large regime of the temperature-independent resistance. The reason for this is likely to 
be due to the important physical effects missed by the model ([T|), namely the presence 
of random induced charge on the superconducting islands. As was shown in a number 
of Josephson junction studies (see e.g. [l^l ) the induced charge on each island exhibits 
very slow random fiuctuations and is therefore inherently random variable. Most likely 
the time dependence of these fiuctuations can be neglected and the induced charge can 
be regarded as a quenched random variable qi that should be added to the Hamiltonian 

The properties of the model ([5]) are not currently well understood; in particular, it 
seems likely (see (l9l |) but was not proven that insulating state has many glassy features 
responsible for intermediate 'normal' phase. It is, however, clear that in this model 
the transition, or a series of transitions occurs at Ej/Eq ~ 1, which corresponds to 
G ^ 1 in the normal state. Near the critical point the gapless excitations correspond 
to collective modes build of electron pairs while the electron spectrum remains fully 
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gapped. In analogy with spin glasses, one expects that effective frustration introduced 
by random charges and magnetic field leads to a large density of low energy states. This 
might explain the observed temperature independent resistance that varies at least by 
one order of magnitude around G ~ 1 as a function of magnetic field. Q 

As we show in section [TT^ both the data and theoretical expectations for models (|ll2p 
differ markedly from the behavior of the homogeneously disordered films. 

A spectacular property of the superconductor-insulator transition in granular mate- 
rials is that a strong magnetic field applied to the system in the insulating regime results 
in a dramatic increase of the conductance. Qualitatively, the reason for this behavior 
in granular systems is that in the absence of the field the single electron excitations are 
absent due to superconducting gap while pairs are localized due to Coulomb energy and 
random induced charges. Large field suppresses superconducting gap, which allows trans- 
port by individual electrons that is characterized by a much larger tunneling amplitude 
and lower (by a factor of four) effective charging energy. This effect was reported by [2^ 
where strong magnetic field was applied to Al grains immersed into Ge insulating matrix 
and giant negative magneto-resistance was observed. Similar behavior was reported for 
homogeneous films of InO deep in the insulating regimepH. This similarity indicates 
that the main reason for this effect, which is that the pairing of the electrons survives 
deep in the insulating state, also holds for homogeneous films of InO. The quantitative 
theory of negative magnetoresistance in granular superconductors was developed in (22j 
for the case of relatively large inter-grain conductances Gij « G ^ 1 , in which case the 
negative magnetoresistance effect is small as 1/G. Recent review of theoretical results 
on normal and superconductive granular systems can be found in [2^ . 

1.1.2. Coulomb suppression of in uniformly disordered thin films. 

The scenario described above assumes that superconductivity remains intact inside 
each grain. An alternative mechanism for the suppression of superconductivity by 



Coulomb repulsion was developed by Finkelstein |24l . |8| , building upon earlier pertur 



bative calculations [25| . Finkelstein effect becomes important for very thin strongly but 
homogeneously disordered films, as well as quasi- ID wires made out of such films [26| . 
Contrary to the Coulomb blockade scenario, the system is supposed to be "uniformly 
disordered" , with no superstructures such as grains coupled together by weak junctions. 
Somewhat similar idea was proposed [27, 28 1 for three-dimensional materials near the lo- 



calization threshold. The essence of Finkelstein effect is that Coulomb repulsion between 
electrons gets enhanced due to slow diffusion of electrons in highly disordered film, which 
results in the negative contribution to the effective Cooper attraction amplitude at small 
energy transfer e: 

A(£) = Ao--^ln-^ (3) 

where g = /i/(2e)^i?n is dimensionless film conductance, Aq is the "bare" Cooper at- 
traction constant defined at the scale of Debye frequency lod, and r* = maxr, t(6//)^, 
where r and / = vpT are the mean scattering time and mean free path, and h is the film 
thickness. The suppression of attraction constant Eq.([3]) leads immediately (we assume 
here ujd ^ 1/'''*) to the result obtained early on in the leading order of the perturbation 
theory ^ 

^^^^-^m^^ (4) 
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The leading terms in the perturbation theory for Tc can be summed by means of the 
renormalization group method developed in [2J] . In the leading order over 1/ g one 
gets 



ln(ft/r,on^T^ 



(5) 



/STrg + ln(fi,/rcoT* 

According to Eq. ([5]), vanishes at the critical conductance gcP = ln^(ft/TcoT*)/(87r), 
(which needs to be large enough for the theory to be self-consistent). At lower (but still 
large compared to unity) conductances, the material never becomes superconducting; it 
stays metallic at least down to very low temperatures Tjoc ~ (/i/r* ) exp(—47r(7) where 
weak localization crosses over into the strong localization [2^ . 

This mechanism of superconductivity suppression, described by Eq.(IS]), might be 
called "fermionic", as opposed to the "bosonic" mechanism discussed in previous sub- 
section [lol |. Within this mechanism, superconductivity is destroyed at relatively large 
conductances QcF ^ 1, thus a direct superconductor- insulator transition does not seem 
to be a natural option. 

The theory of the superconducting-insulator transition outlined above neglects the 
mesoscopic fluctuations of the interaction constant, ^(r). This assumption was questioned 
on phenomenological grounds by Kowal and Ovadyahu [soj . The role of these fluctuations 
become larger when superconductivity is strongly suppressed because as follows from 
in this regime even small mesoscopic fluctuations of gir) lead to a large spatial 
fluctuations of the effective coupling constant A(£,r). In its turn, the fluctuations of 
the effective coupling lead to the local spatial fluctuations of the transition temperature 
Tc{v), which becomes very strong, 5Tc/Tc > 1, for nearly-critical conductance g ~ g^p^ as 
shown in Ref . (sH . Note that mesoscopic fluctuations of Tc remain small in the universal 
case of short range disorder if Coulomb suppression of superconductivity is not taken into 
account, even in the vicinity of the upper critical field Hc2{0) at very low temperatures. 



These results demonstrate the inherent inhomogeneity of superconducting state near 
the critical point at which it is destroyed by disorder. Thus, it seems likely that the regime 
close to the superconductor-insulator transition is described by the effective model of su- 
perconducting islands (appeared due to spatial fluctuations of local attraction constant) 
coupled by weak SNS junctions. The theoretically consistent description of this physics 
is still lacking, the difficulty can be traced to the absence of tunnel barriers separating 
fluctuation-induced superconducting islands from the surrounding media. In the absence 
of these barriers charging effects become non-local which makes the determination of the 
effective Coulomb energy a nontrivial problem, furthermore the presence of a large nor- 
mal part implies dissipative (non-local in time) dynamics of the superconducting phase. 
A toy model of this type was solved in [s^ ; this model describes artificial superconductive 
islands in a good contact with disordered thin film. 

In conclusion, in the fermionic mechanism the Coulomb interaction is enhanced by 
disorder leading to the suppression of superconductivity. The state formed when the 
superconductivity is suppressed is likely to be a poor metal characterized by a large 
resistivity and finite density of states at the Fermi level. 
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1.1.3. Localization versus Superconductivity. 

The third alternative mechanism for the superconductor-insulator transition is due 
to the localization of single electrons. In this scenario the effects of Coulomb interaction 
are neglected, whereas local (in space) Cooper attraction is treated within standard BCS 
approximation. Here we focus on the case of bulk disordered materials or sufficiently thick 
films in which localization remains a three dimensional effect. Abrikosov and Gor'kov [sBj 
and Anderson [s^ have shown that potential disorder does not affect thermodynamic 
properties of usual s-wave superconductors. More precisely, this statement (based upon 
the presence of time- reversal symmetry and called "Anderson theorem") means that 
the parameter Tct/H does not appear in BCS theory as long as magnetic field and/or 
supercurrent are absent. However, localization of single-electron eigenstates at very 
strong disorder leads to appearance of an additional energy scale 5l = l/voL^i where L 
is the single-electron localization length and is the density of states (per single spin 
projection). The meaning of 5l is just average level spacing inside typical volume where 
wavefunction is localized. One expects that at large 5l the superconducting pairing 
between electrons is suppressed. Competition between superconductivity and Anderson 
localization was studied orig inally in mid-80's [i,i,Sl. Their major conclusion was that 
Anderson theorem is valid and superconductivity survives provided that the condition 

Tc » 6l (6) 

is satisfied. The reasoning leading to Eq.® is that for Cooper instability to develop, 
characteristic energy spacing between hybridized Cooper pairs (which are supposed to 
be localized in the same region of size L) should be smaller than typical energy scale 
corresponding to the Cooper instability. On the contrary, no superconducting long-range 
pairing seems possible when level spacing 5l strongly exceeds Tc, in spite of the presence 
of inter-electron attraction (we assume that Cooper attraction survives when the single 
electron states are localized as long as 5l is much smaller than Debye energy ujd)- 

We will show below that the analysis presented in 0, 0, Hi is not complete in two 
important respects. First, the absence of long-range superconductive order does not 
necessarily mean that Cooper pairing is totally negligible; we will show that in the range 
7c ^ <5l ^ Cooper correlations leads to formation of the hard-gap insulator instead 
of usual variable-range one. This gap is of the same origin as the "parity gap" studied 



by Matveev and Larkin 37| in the context of ultra-small grains of good superconductive 
metal. Second, the notion of eigenfunction fractality (that was not known when the 
theoryd i, i] was developed) has to be taken into account and leads to important 
physical consequences. In the present paper we are going to fill both these gaps; it 
will be shown that fractality of electron eigenfunctions changes qualitative features of 
superconductive state, and even modifies the condition Sl ~ T'c for the critical region 
where superconductivity is finally destroyed. 

More recently the issue of competition between localization and superconductivity 
was reconsidered in the important paper by Ghosal, Trivedi and Randeria (Hj. They con- 
sidered two-dimensional lattice model of superconductivity with moderately strong local 
attraction (negative-U Hubbard model) and on-site disorder and studied it numerically 
by two methods: by solving the self-consistent Bogolyubov-de Gennes equations, and 
by solving BCS pairing equations in the basis of exact single-electron eigenstates. They 
demonstrated that with increase of local disorder superconducting state is transformed 
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into the insulating one. The latter possesses sharp gap in the density of states but does 
not show coherence peaks. The energy gap was found to be non-monotonic as function of 
the disorder strength. It was also shown that superconductivity is very inhomogeneous 
in the crossover region, with disorder-generated "islands" of large pairing amplitude. We 
will see below that qualitative features of the results obtained in Rcf. [5| arc very robust 
and survive in a continuum weak-coupling BCS model that we consider in the present 
paper (sec section [2]). The drawback of the treatment developed in Ref. Q is that it does 
not allow the quantitative analysis of the physical properties as function of main pa- 
rameters of the problem (coupling strength A <C 1 and proximity of the Fermi-energy to 
the localization edge, \Ep — Ed <^ Ep), due to limitations imposed by purely numerical 
methods. A major drawback of most conventional numerical methods is their inability 
to study the regime characterized by dramatically different energy scales, in particular 
Tc ^ Ep. Development of the method to study this regime is the goal of the present 
paper. 

The importance of the analytical treatment of the weak coupling regime Tc <C Ep is 
demonstrated in particular, by the numerical work jsif which has studied the 3D disor- 
dered Hubbard model with strong local attraction (4 times larger than bandwidth). In 
this regime the electrons arc strongly bound to each other even in translationally invari- 
ant systems. One expects that the mobility of the formed pairs is less than the mobility 
of the original electrons; this enhances the effect of the disorder. This expectation is 



conformed by the data[38[. We show that in the physically relevant regime of weak at- 
traction, the situation is opposite: the superconductivity survives deep in the regime of 
localized states. 

1.2. Experimental results on S-I transitions. 

Phenomcnologically one should distinguish at least three types of materials that dis- 
play superconductivity suppression with the increase of the disorder: granular systems 
[20I [39I I , nominally homogeneous films that exhibit superconductor- metal- insulator tran- 
sition and homogeneous films that show direct superconductor-insulator transition. j40j 
We shall discuss only the latter class in this paper, the materials that exhibit it are 
thick (more than 20 nm) InO^: films, thinner (abound 4nm) TiN films and very thin (few 
atomic layers) Be films. Recent work reports similar behavior also in disordered epitaxial 
films of NbN with varying disorder [41|. The goal of this section is twofold: to discuss 
the data that allow to exclude Coulomb mechanisms (both fermionic and bosonic) of the 
superconductor-insulator transition in these films and to briefly summarize the data on 
these films that need theoretical explanation. We begin with the first. 

The most direct experimental evidence that allows to exclude the 'fermionic' mecha- 
nism of superconductivity suppression discussed in section [TTT21 in homogeneously disor- 
dered films that exhibit direct superconductor-insulator transition comes from the recent 
tunneling data. In these experiments one observed that the suppression of the supercon- 
ductivity either by disorder or temperature is not accompanied by the su ppr ession of the 
gap, which remains intact or even increases reaching 2Ai/Tc G 6 — 9. [43 Instead, as 
the insulator is approached, one observes the disappearance of the local coherence peaks 
which, in addition, vary strongly from one point to another ( see also Ref. [4^) . This 
behavior (and especially the temperature dependence of tunnelling conductance) is in a 
striking contrast to what is expected for the fermionic mechanism. Less direct evidence 
is provided by the data 3, 45 1 showing that superconductivity exists up to a very strong 
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disorder corresponding to g ^ 1, which is at least a factor of two smaller than the one 
expected in the fermionic mechanism ([5|). 

As explained in section 11.1.21 the basis of the fermionic mechanism is the idea that 
Coulomb repulsion is enhanced by disorder which results in the effective suppression of the 
attractive interaction that leads to superconductivity. The actual equations arc derived 
in the assumption that bare Coulomb repulsion is very strong but is reduced by screening 
to the universal limit in which the effective Coulomb repulsion constant is equal to unity. 
The dimensionless parameter characterizing the strength of the Coulomb interaction and 
its screening is 2a/(Tcii), where a = {e^kp /Q-k'^) {kpl) is the residual conductivity and k 
is the dielectric constant due to electrons far from the Fermi energy {\E — Ep\ > up)). 
Coulomb interaction is effectively strong provided that 2a/{TcK) ^ (Co/ascr)^/^ S> 1, 
where Cgcr is the Thomas-Fermi screening length, is the coherence length in a dirty 
superconductor. If instead 2a/{TcK) <i; 1 the coefficient in front of the logarithm in 
piSp becomes small so that enhancement of Coulomb repulsion become important only 
at exponentially low energy scales. In a very dirty metal (such InOj^film) with a short 
mean free path kpl '--^ 0.3 and low carrier density e^kp ^ 5000K (see Ref . j46l| ) the ratio 
tj/Tc ~ 10. Thus, the effects of Coulomb interaction in these films become unimportant 
if dielectric constant k ^ 10. The direct measurements of the dielectric constant deep 



in the insulator regime give k > 30 j47l |. one expects that it can be only larger in the 



vicinity of superconductor- insulator transition, so 2a/{TcK) <^ 1 in these films which 
makes the fermionic mechanism irrelevant. 

The microscopic origin of this large dielectric constant is likely to be due to the low 
density of the carriers, ^ lO^^cm^"^ in these conductors, and the peculiar structure 
of their density of states in which the Fermi level is located in a large dip.ji^l. In this 
situation the density of the electrons distant from the Fermi level, \E — Ep\ > ujp), is 
high, which results in a large screening of the Coulomb interaction. 

Now we turn to the possibility of the Coulomb driven transition similar to the one 
of Josephson arrays fsection ll.l.ip . This physics is due to the long range nature of the 
Coulomb interaction; the estimates of the associated energy scales given below show 
that Coulomb interaction is sufficiently well screened even in poor conductors so that 
the corresponding energy is too small compared with all other energy scales; this rules 
out this possibility. These arguments are quite general and apply to other effects that 
originate from the long range part of the Coulomb interaction. 

We begin with energy scales that are known to be relevant experimentally for this 
problem. Namely, the superconducting gap on the ordered side and activation energy 
on the insulating side. These energies are large: the superconducting gap is around 
A ^ hK while the activation gap is even larger Tq ^ \Q — IbK (see below). We now 
estimate the Coulomb interaction at the scale of the superconducting coherence length, 
^0- Because Tc/Ep < 10"'^ in these materials, this length, even for a very poor conductor, 
cannot be too short: ^ IOtito. Large value of the dielectric constant in the parent 
insulating compounds, k > 30 , implies that in the absence of the mobile electrons with 
\E — Ep\ < up), the effective charging energy at scales would be Eco = e^/^^o ^ 50 K. 
Screening by conduction electrons with energies Tc < \E — Ep\ < ujd decreases it further. 
To estimate this effect we note that at scales less than the properties of the electrons 
are similar to those at the mobility threshold. At the threshold, the dielectric constant 
grows with scale, L, according to the scaling law (LjlY with x > 1. This results in the 
additional factor (VCo)^ in the effective Coulomb energy at scales which reduces it to 
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Ec < IK. Thus, the effective Coulomb energy is much smaller that all relevant energy 
scales and cannot be the driving force of the transition. Note that this estimate becomes 
incorrect in the presence of thin insulating barriers between the grains, thus allowing 
for Coulomb driven transition in inhomogcneous materials. The absence of structural 
inhomogencitics that might lead to such barriers in films of InOa; was shown in (30| ; later 
studies[49| also reported the absence of inhomogeneities in TiN films. 

Another argument against Coulomb driven transition is provided by a completely 
different phenomenology of the transition in the films and in Josephson arrays: in the 
former one observes direct transition to the insulating state characterized by a large gap 
and activation behavior of resistivity, the transition can be driven either by the increase 
of disorder or by magnetic field. In contrast, in Josephson arrays the transition driven 
by the field is characterized by a large intermediate regime of temperature independent 
resistivitv. (ol. [Hoj . Furthermore, there is no reason to expect the disappearance of cohcr- 
ence peaks in some places and not in others as one approaches the transition (reported 
in (42|) in the array of superconducting grains. 

Finally, Coulomb mechanism cannot explain the large value of the activation energy 
in the insulating state in the vicinity of the transitions. One expects that when charging 
energy becomes large enough to suppress transport by Cooper pairs, single electrons 
should dominate transport which implies activation energy equal to superconducting 
gap. Instead one often observes large activation energies: Tq « ISX in InO^; films (see 
Fig.l of (llj) and even larger in Be ones lsil. 

We now summarize (see also the review [40[) the phenomenology of the direct superconductor- 
insulator transition in homogeneous films as exhibited by three different systems: thick 



amorphous InO.^ films [4i,[30l 21, 44 [5l[53|, thin TiN films [Sj, |45|, |55| and extra-thin 
(below 1 nm) Bcrillium films 56|, ISll |57| ; somewhat similar phenomena were observed 
recently in the patterned Bi film with honeycomb array of holes [5^ . 



1. 



2. 



3. 



5. 



On insulating side of SIT, low-temperature resistivity curves show simple acti- 
vated behavior, R(T) oc exp(Tn/T), which crosses over into Mott [s^], R(T) cx 
exp(TM/T)i/4, or Efros-Shklovsky ©[Hll, R{T) oc exp{TEs/Ty^^ variable-range 
hopping at higher temperatures. This behavior is highly unusual: in hopping insu- 
lators where the activation is frequently observed at high temperatures it crosses 
over to some fractional (variable-range) behavior upon the temperature decrease. 

At high magnetic fields the films on both sides of the SIT show large negative 
magnetoresistance 2l|, |4J, |52|, |53|, |57|, |55|, |45| . It is important that such behavior 
was observed even in films characterized by a very high activation energy. 

At low fields all insulating films close to SIT show positive magneto-resistance at 
low fields [H [Hi Si [H [Sg. 

Ultra-low-temperature measurements on nearly-critical samples of a-InO^; [H^l and 
TiN (45! revealed a very sharp jump (by several orders of magnitude in current) in 
nonlinear I{V) curves. 

The resistance of Be [s^ and TiN films [H^ approach the quantum resistance h/e^ 
at very strong magnetic fields and low temperatures. Very recently the importance 
of Zeeman pair-breaking for the properties in this regime was demonstrated in 
Ref. [13 . 

The properties of the quantum critical point that separates superconductor and in- 
sulator are not fully established. The samples corresponding to the critical disorder 
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which separates sup erconducting and insulating behavior display insulating behav- 
ior of the i?(T) [61, 45, 5l|; this suggests that the disorder driven superconductor- 
insulator transition is not described by a self-dual theory proposed in [l^ or that 
the critical regime where this behavior sets in is very narrow. Weakly supercon- 
ducting samples can be driven into insulating state by the application of magnetic 
field. Scaling (or a lack of thereof) near this quantum critical point is subject 
of controversy in both the value of the critical resistance and the scaling expo- 



nents. Some works[55|, |61|, |42| report critical value of the resistance larger than 
Rq = 6.5kD, and exponents that do not agree with the theoretical predictions [l7j. 
In contrast, recent paper 6^ reports both the critical value of the resistance and 
exponent in a perfect agreement with the theoretical predictions based on dirty 
boson scenario [17|. 

These results can be summarized by the low temperature phase diagram in {H, G) 
plane sketched in Figure [T] 

It is tempting to explain these observations by the scenario in which the disorder 
destroys the global superconducting coherence while preserving the local superconducting 
gap, A, even in the insulating phase. Indeed, in the absence of pair coherence the 
conductivity is due to the thermally excited fermionic quasiparticles. The density of 
these excitations would be additionally suppressed by local superconductive gap A that 
is formed at low temperatures: ni{T) ^ exp(— A/T) . This would explain the crossover 
to activation behavior at low T. High magnetic field suppresses local A{H) thus leading 
to very large negative magnctorcsistance. Low-field positive magneto-resistance could 
be then associated with frustration induced by magnetic field which eliminates the last 
vestiges of superconducting coherence thereby shifting the array further into insulating 
side (whereas not yet affecting local superconducting gap). 

This scenario would be naturally realized in the granular material where supercon- 
ductivity remains intact in each grain whereas global coherence appears only due to 
Josephson coupling which competes with charging energy (see section [T.l.ip [63[. How- 
ever, this assumption of hidden granularity of the SIT films is not plausible for many 
reasons explained above. 

The correct theory should be also able to explain the large value of the activation 
energy on the insulating side of the transition. Assumption that it is due to the local su- 
perconducting gap in weakly coupled grains is not sufficient because this energy is higher 
than the gap in a less disordered film of the same series that shows superconductivity. 
Indeed, maximal observed values of Tq were about 15K in work[2lj and lli^ in workjsO] 
which are significantly larger than maximal superconducting gap A « 5K shown in Fig.l 
of Ref. (note that samples studied in [2l| were from the same source as in iil. SOj). 



Furthermore, this assumption would lead to the conclusion that the transition is due to 
competition between Coulomb and Josephson energies which are both much smaller than 
To ; because Ej = GA/2 this would imply that G <C 1 in a direct contradiction with the 
data. 

Large value of the activation energy, together with inefficiency of electron-phonon 
coupling at low T provides a natural explanation for a large jump in I{V) observed 
in the insulating state f63|. T he detailed predictions of the theory [g^] were recently 
verified experimentally in |65| . The important ingredient of this theory is a hard gap 



for single electron excitations, similar to the one directly observed in 42[ and inferred 
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Figure 1: (Color online) Sketch of the experimental phase diagram of homogeneously disordered films 
(upper panel) at T — >■ and its interpretation in the theory of superconductivity developed in the 
electron system close to mobility edge (lower panel). 
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from the resistivity data. Thus, the phenomena of I{V) jumps does not impose additional 
constraints on microscopic theory. Moreover, very similar jumps were observed previously 
in the system which seems to have nothing to do with superconductivity [66| . 



Real challenge to the theory is presented by the magnetoresistance data from Ref. [21| : 
straightforward interpretation of the negative magnetoresistance data as being due to 
superconducting gap suppression in individual grains leads to unphysically large values 
of critical field needed to destroy the superconductivity in such grains. For instance, 
the field of 8 T was observed to produced only moderate ( R{H = 8T)/i?(0) sa 0.5 ) 
negative magneto- resistance in a sample characterized by To w 15 K (as determined in the 
temperature range 1.3—5 K). Interpreting this effect as being due to the suppression of the 
pairing gap, we find Tq — Tq{H = 8T) « 0.7K, which is about 5% of Tq. Interpolating this 
dependence we find that the field necessary to destroy completely the superconductivity 



in each grain is huge: H^^^ '--^ 50 — 80 T 



iry 

i 

Such large values of the critical fields are impossible for realistic grain sizes. Indeed, 
critical orbital magnetic field for a small (radius R < where ^ = y/HD/A is the 



coherence length) superconducting grain is [671 122| 



rrest ^ ^OOT 

where R and ^ are measured in nanometers. Using a typical diffusion constant for a poor 
metal D « Icm^/s and allowing for a very high gap value A = 10 K, we find ^ = 8.5 nm. 
Together with the lowest bound for the grain radius R = 6 nm in which the distance 
between the levels does not exceed the superconducting gap 

6 = {AiyoR^y^ < A 

it leads to H^g* < 20 Tcsla which is still much smaller than H^g"^ above. These estimates 
did not take into account the spin effect of magnetic field that would further decrease 
the value ofmf. 

We conclude that the observed activation gap cannot be explained as BCS gap in small 
grains composing the material: it is too wide and too stable with respect to magnetic 
field. 

A number of works argued that mesoscopic fluctuations might lead to the appear- 
ance of inhomogeneous superconductivity (self-induced granularity) in the vicinity of the 
transition even in the absence of structural granularity (sol . |40| . |45| . The direct com- 
putation shows that these speculations are correct for the fermionic mechanism of the 
superconductivity suppression in two dimensions [sij . We expect that the self- induced 
granularity that appears due to this mechanism does not lead to thin insulating barriers. 
It is therefore characterized by a small value of the Coulomb interaction between the 
'grains'. Thus, it can be ruled out as a mechanism of direct superconductor-insulator 
transition in homogeneous films by the energy scale arguments given above. 

An important unresolved issue is the nature of carriers responsible for the transport 
in the insulating state: are they Cooper pairs or single electrons? One expects that the 
presence of the superconducting gap in the insulating state implies that the transport is 



Data of Ref. [s^ show that 32 Tesla field is not sufficient to fully suppress negative magnetoresistance. 
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dominated by Cooper pairs in the vicinity of the transition and this was indeed observed 
in ultrathin Bi films |68l|. However, one expects that the transport is dominated by 
single electrons further in the insulating state where activation behavior was observed. 
Unfortunately there are no data to confirm this. 

To summarize: experimental data on SIT in amorphous materials call for a new 
mechanism of a gap formation, which is somehow related to the superconductivity, but 
is different from the usual BCS gap formation. In the vicinity of SIT this mechanism 
should lead to a "pseudo-gap" features in R{T) behavior and tunneling data. 



1.3. Main features of the fractal pseudospin scenario. 

The fractal pseudospin mechanism (briefiy presented in j69| ) should be viewed as 
alternative to both 'boson' and 'fermion' scenarios. Here we argue that it is fully com- 
patible with the data. The key elements of this approach are in fact quite old: (i) An- 
derson's reformulation [t^ of the BCS theory in terms of "pseudo-spins", (ii) Matveev 



Larkin theory [37| of parity gap in ultra-small superconducting grains and (iii) frac- 



tal properties of single-electron eigenfunctions with energies near the Anderson mobility 
edge [2l|,[z2,[2i- 

Qualitatively, in this scenario the electrons near the mobility edge form strongly 
coupled but localized Cooper pairs (notion first introduced in [2l|, |6lj, see also [^) due 
to the attraction of two electrons occupying the same localized orbital state. These pairs 
are characterized by a large binding energy which is responsible for the single electron gap 
To observed in transport measurement in the insulating state. At temperatures below Tq 
the system can be described as a collection of Anderson's S = pseudo-spins, whose Sj 
components measure the Cooper pair occupation number and components correspond 
to pair creation/annihilation operators. Superconductivity in this system is due to the 
tunneling of Cooper pairs from one state to another. It is essential that it competes 
not with the Coulomb repulsion but with the random energy of the pair on each orbital 
state. In spin language it is described as a formation of non-zero averages (S^) due 
to "off-diagonal" S~S^ + S^SJ coupling in the effective Hamiltonian which competes 
with random field in z— direction term hiSf . Large values of the binding energy and off- 
diagonal interactions are due to the properties of localized nearly-critical wave functions; 
the main features of these wavefunctions are their strong correlations both in real space 
and in energy space, and their sparsity in real space (see section [^T^ . The resulting 
phase diagram is shown in Figure [TJ 

The fractal pseudospin scenario has many common features with the bosonic mecha- 
nism, but it is distinct from it in a few important respects: 

• pseudo-gap energy scale Tq is independent from the collective energy gap A 

• fractal nature of individual eigenstates implies a large " coordination number" Z S> 
1 of interacting pseudospins away from the superconductor-insulator transition. 
Close to the transition Z drops, resulting in very inhomogeneous superconductive 
state and an abrupt decrease of Tc. 



distribution of superconducting order parameter in real space is extremely inho- 
mogeneous, thus usual notion of space- averaged order parameter A is useless even 
qualitatively, and the "Anderson theorem" is not applicable. 
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In the main part of the paper we present theoretical arguments in support of this new 
scenario. We restrict our discussion to the three-dimensional problem which is appropri- 
ate for the electron wave function behavior in most films. It is possible that the physics in 
the near vicinity of the transition is dominated by large scales where the two dimensional 
nature of the films become important, the details of the crossover to this critical regime 
is beyond the developed theory. We will assume below that localization effects are not 
very strong, allowing for the presence of phonon-induccd attraction between electrons. 
Clearly, the necessary condition for that is 6l ^ wd- 

The theory that we develop starts with the single electron states of the non-interacting 
problem, so it is not applicable to describe the physics in high magnetic fields where these 
states change significantly. Thus, interesting physics of the metallic state with resistance 
approaching h/e'^ is beyond the applicability limits of our theory. 

2. Model 

2.1. BCS Hamiltonian for electrons in localized eigenstates. 

Wc consider simplest model of space-local BCS-type electron-electron attraction, 
Vint = 5<5(r). It is assumed, as usual, that this attraction is present for electrons with 
energies E in the relatively narrow stripe E Cz Ep ± uju around Fermi energy. However, 
we will see below that in contrast with the usual BCS theory, the parameter ljd will not 
enter our final results. The Hamiltonian represented in the basis of exact single-electron 
eigenstates ipj{r) becomes 

jcr i,j.k.l 

where 

A%,,- J drtl;*{r)rj{r)Mr)Mr), (8) 

= Ej — Ep is the single-particle energy of the eigenstate j counted from Fermi level, 
Cja is the corresponding electron annihilation operator for the spin projection cr, vq is the 
density of states (per single spin projection) and A = gvo ^ 1 is dimcnsionless Cooper 
coupling constant. 

Note that writing Hamiltonian in the form ([7]) wc omitted the Hatree-type terms 
which do not contribute directly to the Cooper instability. Such terms are known to be 
negligible when single-electron states are extended (see discussion in Ref. 0]); the issue 
of their importance for critical and, especially, localized states is more delicate. Below in 
section 16.11 we present results for the superconducting transition temperature obtained 
with and without account of the Hatree-type terms. The comparison, shown in Fig. 
demonstrates that these terms, while changing the quantitative results somewhat, do 
not affect our main qualitative conclusions. Therefore, in order to keep the arguments 
as simple as possible, we neglect Hatree terms in the main part of the following text, d 



^However, those terms are important for a quantitative description of superconductivity in the region 
of locahzed single-particle states, as shown in Ref. where 2D problem in the limit of strong disorder and 
strong attraction was studied numerically. In particular, they might lead to additional inhomogeneous 
broadening of the coherence peaks observed in tunneling experiments, see section r6.2.3l 
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Unless specified, we will not consider magnetic field effects, thus eigenfunctions ipj{r) 
can be chosen real. In the following we will use frequently a simplified Hamiltonian ([7]) 
where only pair- wise terms i = j and k ~ I are taken into account: 



where 

Af,fe = 1 dvi,]{v)^l(v) . (10) 

Eq. is the minimal Hamiltonian that includes hopping of pairs necessary to establish 
a global superconducting order. It plays the same role for our theory as the BCS Hamil- 
tonian with « = p, J = — p, fc = p', Z = — p' for usual theory of superconductivity. We 
will discuss the accuracy of this approximation below in Sec. 14.21 

2.1.1. Ultra-small metallic grain. 

Here we rederive the known results for the model (|9]) applied to ultra-small metal 
grains; this derivation will provide the starting point for our solution of the model Hamil- 
tonian ([7]) or ([5]) for the electrons with Fermi level near mobility edge. 

Pairing correlations in metallic grains of very small volume V, with level spacing 
S ~ (vqV)^^ comparable to the bulk superconductive gap A were considered in many 
papers, see review 1741 . The issue which is most important for this work is the parity 
gap introduced in [37| to characterize pairing effects in ultra-small grains with (5 ^ A. 
The work [s^] assumed the simplified Hamiltonian © with identical matrix elements 
Mjk = 1/V. 

More complete treatment of a weak electron-electron interaction in small metallic 
grains is given by [ts! ] where it was argued that to the leading order in the small pa- 
rameter tqS, where tq is the flight (or diffusion) time for electron motion inside grain, 
all off-diagonal terms Mijki can be neglected. In the relevant terms the indices of Mijki 
should be pairwise equal. Because for small grains the wavefunctions tpj (r) are essentially 
random Gaussian variables subject only to orthogonality and normalization conditions, 
the matrix elements pop appearing in the Hamiltonian ^ are given by 

Mj^k = ^ Af,,^M, = |, (11) 

so that the full Hamiltonian can be expressed in terms of the total number of electrons n, 
the total spin S and the operator T ~ CkiCk-\ related to Cooper pairing correlations: 

- XST^T (12) 

It is essential for the validity of ([T^ that all matrix elements Mjk with i j are equal 
to 1/V, while diagonal terms are three times larger, only in this case it is possible to 
represent Ea. (jl2p in terms of the total density, spin and pairing operators. 

Attractive interaction implies that S* = in the ground state. Because n is conserved 
for isolated grain the properties of this model in n = 0, S ~ sector are equivalent to 
the properties of the simplified Matveev-Larkin model ^ which takes into account only 
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the last term in ([T2|). The first term in ([12)) is important for the correct evaluation of 
the coefficient of the interaction term with j = fc in Eq.Q because only 1/3 of it should 
be assigned to the interaction in the Cooper channel, since other 2/3 contribute to the 
" n" and " S" terms of the Hamiltonian. 

In the limit (5 3> A one can use the perturbation theory with respect to pairing 
Hamiltonian (fT^ . In the lowest order in A, neglecting all terms except diagonal ones, 
one finds that the energy of two identical grains with even number of electrons, n ~ 2k, 
and zero spin is by AE = 3X5 less than the energy of the same two grains with 2fc+ 1 and 
2k — I electrons and spin 1/2. [f| It means [s^l that the average ground-state energy of a 
grain with even number of electrons is lower by "parity gap" Ap = |A(5 than the energy 
of the same grain with odd number of electrons. Note that Cooper pairing contributes 
1/3 of this energy difference. This result is valid only in the limit of a very small coupling 
constant A, when all the terms with j 7^ fc in Eq.® can be neglected. In a more general 



case these terms must be taken into account which leads 371 to the renormalization of 



the coefficient of the T^T term in the Hamiltonian (|12p which becomes 

XR = X/{l-X\n{ojD/S)). (13) 

After this renormalization the coefficient of the Cooper pairing becomes dominant. In- 
troducing bulk energy gap A = ujoe~^^^, one finds [37| parity gap which is valid for all 
A<5: 

A, = j^ + A (14) 

where the second term is due to the first term in the Hamiltonian (jl2p and is small 
compared to the main term. The result (jl4p shows that parity gap grows with the 
decrease in the grain size. Note that parity gap would not appear if one does 
note take into account the double-diagonal terms in ([3]), which are totally irrelevant in 
the usual BCS theory of bulk superconductivity. Below we will find somewhat similar 
behavior in the case of bulk Anderson insulators. 

2.1.2. Vicinity of the mobility edge. 

In the bulk Anderson insulator with the Fermi energy near the mobility edge, the 
typical energy scale replacing S is 

5l = l/ii^oLU, (15) 

where Lioc is the localization length. It was argued in Q that localization is irrelevant for 
superconductivity if 3> In the opposite limit Tc <C (5l pairing correlations between 
electrons localized on different orbitals are irrelevant. Localization length depends on the 
Fermi-energy (in the scaling region Lioc 3> £) as 

(16) 



^For different grains containing different number of particles one needs to take into account different 
cfiemical potentials in these grains but the final conclusion remains unchanged. 
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Figure 2: Inverse localization length as function of proximity to the mobility edge, obtained numerically 
for 3D Anderson model with Gaussian disorder of width W = A. The values shown here were extracted 
from the numerical computation of the inverse participation ratio in this model and its conversion into 
the localization length by Eas. lll5l71ll . Full line is a fit to l/Ljoc = 1-87 ■ (E/Ec ~ l)^'^. 



where Ec is the position of the mobihty edge, v is the locahzation length exponent. Nu- 
merical data [7^ show that in a very narrow vicinity of the mobility edge {{Ec — Ep)/ E^ <^ 
1) of 3D Anderson model the localization length obeys the scaling dependence with 
V ~ 1.57. As we shall see below, the range of energies relevant for the superconductor- 
insulator transition is relatively wide, {Ep — Ec)/E^ < 0.5; in this broader range the 
localization length follows the same scaling behavior (|T6)) but with as somewhat differ- 
ent exponent f « 1.2. We illustrate this by Fig. [5] that shows the localization length 
obtained for 3D Anderson model with Gaussian disorder (see section 12.21 and Eq. 
below). The parameter £ in Eq. p^ is the short-scale cutoff of the order of the elastic 
scattering length. The associated energy scale 

Eo = l/(j^o^') (17) 

depends on the microscopic details of the model of disorder and can be small compared 
to Fermi-energy Ep (see next subsection for the discussion of this issue). 

We will assume that Fermi energy is not too close to the mobility edge so that 

Er — Ep Tr 

t^^—^»-^ (18) 

The condition (fT8)) means that localization properties of eigenstates ipj (r) do not vary 
appreciably within the energy stripe Ep ± mainly responsible for development of 
superconducting correlations. Violation of the condition (|18p implies the absence of the 
particle-hole symmetry in the superconducting state; we expect anomalous Hall effect in 
the superconducting state and near transition in this regime. Using Ea. (|16p we find 

5l = Eot""" (19) 



Because 3i^ « 4 ^ 1, the condition ([18)) is compatible with t <C 1 in a wide range of 
parameters which include both small and large ratios 5l/Tc. In the limiting case 5l ^ 
the scale set by localization is much larger than the one set by superconductivity, so all 
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relevant statistical properties of matrix elements pOj) can be computed at i = 0. We 
will refer to this case as the critical regime. We will see in section |4] that deep in the 
limit 5l/Tc 0, the transition temperature approaches its limiting value, which we 
denote as r°. In contrast with the conclusions of Ref. we will find that may differ 
substantially from the usual BCS value T^'-'^ = u}r)e~^^^ for the metal with the same 
value of the Cooper attraction constant. Moreover, we find that the values of T° are 
typically larger than T]^'^^ for weak couplings A <C 1. This unexpected result is related 
to the fractality of electron wavefunctions with energies close to the mobility edge. 

2.2. Fractality and correlations of the wave functions near the mobility edge. 

The exact single-particle eigenfunctions ipj (r) and eigenvalues Ej that enter the model 
Hamiltonian Eq. ©should be found from the single-particle Hamiltonian with disorder. 
The conventional models of disorder are the continuous model of free electrons in a 
Gaussian random potential U{r): 

Hia = ^^+Uir), (20) 

or the tight-binding model with on-site energies £„ being random variables with the 
probability distribution Vden}) — IlnPC^n)- "^^^ latter is known as the Anderson 
model, it is described by the Hamiltonian 

iJib = ^ £„ a|^a„ - ^ aJifln+a- (21) 

n n.m— n+a 

The most common choice of the distribution function p{x) are the box distribution 

^^^^"1 0, if x> \W/2\. 



or the Gaussian distribution 



p{e) = — exp 



/2nW 



2W^ 



(23) 



Increasing the disorder parameter W in the 3d Anderson model pTj) . ((22|) at a fixed 
Fermi energy Ep results in the Anderson localization transition at the critical disorder 
W — Wc (at Ep = the critical value of disorder is Wc = 16.5 for the box distribution. 
Eg. [^^ . Alternatively, the localization transition occurs when Ep is increased at a fixed 
disorder W < Wc beyond the mobility edge Ec] in the following we will mainly use the 
Gaussian model, Eq. ((23)) . The same type of transition takes place in the continuous 
model Eg. ipnj) . The changes in the statistics of wavefunctions resulting from this tran- 
sition do not merely reduce to their localization. Well before all wavefunctions become 
localized they acquire a certain structure where a wavefunction occupies not all available 
space but a certain fractal inside the correlation radius icon - The global picture of an 
extended wavefunction resembles a " mosaic" made of such pieces of fractal with the char- 
acteristics size icon- This peculiar phase (called the "multifractal metal" in Ref. fll\ ) 
appears in the vicinity of the Anderson transition. It persists down to relatively weak 
disorder as long as the decreasing correlation length Lcorr exceeds a microscopic length £ 
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which has a meaning of the minimal length (a pixel) of the fractal structure. In Anderson 

1/3 

model with the box probability distribution the length i « aWc > 2.5a, where a is 
the lattice constant; fractal effects disappear in this model at < 3 ^ Wc only, see 
Ref. [t^]- For the continuous model defined by Eg. ipH)) it is of the order of the elastic 
scattering mean free path. 

As one approaches the mobility edge or the critical value of disorder, the correlation 
radius Lcorr diverges so that the critical wavefunctions are pure fractal (or, strictly speak- 
ing multifractal [77|). On the localized side of the transition the wavefunctions inside 
the localization radius L\oc resemble the one inside an element of the mosaic structure 
of the multifractal metal. This "multifractal insulator" 77[ exists in the vicinity of the 
Anderson transition and becomes an ordinary insulator at strong disorder when Lioc < ^■ 



2.2.1. Wavefunction correlations at the mobility edge: algebra of multi-fractal states. 

We start by describing the multi-fractal correlations of wavefunctions exactly at the 
mobility edge. To avoid confusion we note that for a finite 3d sample of the size Lx Lx L 
the mobility edge is smeared out. The critical multi- fractal states live in a spectral 
window around Ec of the width SE oc L~^/'^, where v is the exponent of the localization 
(correlation) length Lioc(icoir) oc \E — Ec\~'^, such that the value of Lioc(-Z^corr) inside 
this window is larger than L. The number of single-particle states in this window is 
proportional to L^'-^~3iy\ Because v is definitely larger than -i (in fact, Harris criterion 
tells that ly > 2/d — 2/3) it tends to infinity as L — >■ oo, . 

There is a vast numerical and analytical evidence [78[ that the critical wavefunctions 
at the mobility edge obeys the multifractal statistics. This can be seen, for instance, in 
the behavior of the moments of the inverse participation ratio: 

The moments describe the effective volume occupied by the the wave function. At 
the mobility edge they scale with the size of the sample 

{Pq) - r^d-d,)(q-l) ^-d,{q-l) ^ (25) 

where c?g < 3 is the corresponding fractal dimension. For the 3d Anderson model of the 
orthogonal symmetry class (real Hamiltonian) we obtain by numerical diagonalization 
the following results for the first two fractal dimensions: 



d2« 1.29 ±0.1, 0.72 ±0.1. (26) 

The fact that the fractal dimensions dq depend on the order of the moment q implies 
the multiractality of the wave functions. The scaling arguments show that such behavior 



*The fractal dimensions for bigger sample sizes have been studied recently by Rodriguez, Vasquez and 
Roemerlli. They have found d2 = 1.24 ± 0.07, di = 0.63 ± 0.07, 4^"" = 1.35 ± 0.07, ^ = 1.02 ± 0.2. 
They also point out on a large systematic error for et'^'^ related with the finite-size effect. In view of 
the fact that the critical Qc ~ 2.1... 2. 2 is close to 2, the typical d!'^'^ should be found from the condition 
Eq.(|33]l. This gives an estimate d^l^ = 0.84 ± 0.04. 
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of {Pq) implies the power-law correlations of wavefunction amplitudes at different space 
points: 

C,(0,r)^(|V',(r)p^|^,(0)p^)~L-2^''(W' (^^) "\ (27) 
where £ < r < L, and the exponents arc equal to 

ag = d2gi2q^l)-2dg{q-l), (28) 
= 2{q-l)id-dq) (29) 

Note that the sign of aq is positive provided that the moments Pq are only moderately 
fluctuating, so that the scaling behavior of (P^) and (Pq)'^ is the same. This follows from 
the inequality 

and the definition of the fractal dimensions Ea. ipS)) . However, from Ea. ipS)) is follows 
that 

a2 = 3<i4 - 2<i2 « -0.43 ± 0.5. (30) 

Although the error bars are rather large, it is likely that a2 is negative. 

This means that the second moment P2 is strongly, not moderately fluctuating, 
and the L-scaling of (P^) is different from that of (-Pj)^ hi agreement with the early 
coniecturefsoj. Indeed, our numerical simulations on the 3D Anderson model of the 
orthogonal symmetry class show that 

{Pi) OC L-2.16±0.1^ ^p^^2 ^ ^-2.58±0.1_ 

This is consistent with the observation [8l[ that the distribution function of the second 
moment V2{P2) oc P^^^^ ^ power-law tail with p2 ~ 2.6 at relatively large values of 
P2 ^ (P2) ^ L~'^^ (which is cut at P2 = 1). As a result the average (Pf) is dominated 
by the far tail where the one-parameter scaling 7^2(^2) = /2(^2/(^2)) is no longer true. 
The distribution function 7^4 (P4) of the fourth moment P4 possesses even stronger tail 
with V4 OC P4"^ ° at P4 > (P4). In this case the average (P4) is considerably contributed 
by the rare events in the far tail of the distribution. 

In the situation where the assumption of moderate fluctuations of Pq no longer holds 
and the rare events are important for averaging [ts} the correct physical quantity is 
typical average (P5)typ = exp[(lnPg)] instead of the usual one and the corresponding 
fractal dimensions should be defined by 

dT (9 - 1) = -dln{Pq)typ/dlnL. (31) 

Unlike the usual average, the typical average is determined by the body of the corre- 
sponding distribution function Vq{Pq) provided that the tail exponent pq is larger than 
1. Computing the typical averages (Pg)typ and using Ea. ipT]) we obtained (sec also a 
footnote ^ and Ref. [TOj for comparison) for the orthogonal 3D Anderson model: 

4^P K. 1.40 ±0.1, sa 1.08 ±0.1, (32) 

a^yP = 3<P-24''P« +0.44 ±0.5. 
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It is interesting to note [78| that if the tail exponent for a given q-th moment lies in the 
region 1 < pq < 2 (i.e. the typical fractal dimension is determined by the body of the 
distribution for the corresponding moment Pq but the averaged moment is dominated by 
the rare events) the dependence of the typical fractal dimension on q is linear, [f] 

dlyP{q-l) = qa+, q > q,, (33) 

where is a constant. 

Substituting Ea. ([55| into Eg. ([25)1 we immediately obtain: 

afP = (2g-l)4^/-2(g-l)d*/P = 0, q > q,. (34) 

The equality is valid if the condition q > qc is satisfied for both terms entering (pS)) ; 
in all other cases a^^P > 0. Given also that for g > 1 d2q < dq < d, we obtain: 

d > P > 0. (35) 

For 3D Anderson model of orthogonal symmetry class qc ~ 2.2 ±0.1 is very close to (but 
larger than) 2 (see Ref. [?§]). This means that in a typical sample the r-dependencc of 
the correlation function Eg. ipTj) is such that the integral over r is dominated by large 
distances. This will be important for the estimations of the matrix elements done below. 

Long expressions Eas. (|28l29p for the exponents aq and (3q reflect few simple rules: 
Rule(i): for r ^ L the two wavefunctions are statistically independent and the result 
can be foimd using Eq. 



(IV^j(O)l^')^ [L-'^{Pq)f - i-2d,(9-l)-2d 

Rule (ii): for r ^ £ the result should be the same as for the coincident points r = 0, i.e. 

(|V,(0)|'^9) = L-\P2q) ^ L~d,,{2q-l)-d_ 

These rules can be generalized to the case of correlation of different eigenfunctions 

Cqiu^,r)^mr)\^^\^,m'^) 

with the energy difference Ei ~ Ej = uj. In this case the eigenfunctions become statisti- 
cally independent for |r| > i^, where 

Lt - — . (36) 

With this correction, applying the rules (i) and (ii) one can show that at a; > (5 the 
only change in Ea. ((27)l is that the |r|-dependent factor becomes (L(j/|r|)'^~"'' instead of 
(L/|r|)''-"<'. 



^This follows from the fact that the spectrum of fractal dimensions given by the Legendre transfor- 
mation /{«) = —Tq + qa, a = , 'Tq = (l ~ l)'^? determines the scaling with the system size L of 
the number of sites M oc L-'''") where |-(/)(ri)p oc L^" . As a typical event cannot occur at a number of 
sites M < 1 (which is only possible for a rare event), found from d^^^ must be either positive (if 

{Pq) = (P(j)typ) or zero (if (Pq) is dominated by the rare events and is different from {Pq)typ)- In the 
latter case the Legendre transformation implies [1 — gd/dq]{(q — l)d^^^} = 0, leading to the fact that 
(g — l)d^q^ is proportional to q. 
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This allows us to estimate the local averages of different wave functions. In these 
averages the power law correlations will be seen in the energy space. Indeed, substituting 
|r| ^ £ in {L^/\r\)'^~°''' one obtains: 

C,{uj,0) = i.^'Y.(\M0)\''\i^jm"^S{E-E,)6{E'-E,)) 

~ L~''"^{L/ir^ (^^y\ 1,-1-^- (37) 

where {voL''')~^ = S < uj < Eq = {i^q£'^)^^ . For uj < S the correlation function Cq{uj) 
saturates. 

Here, essentially the same rules work: 

(i) at large energy separations uj ^ Eq the wavcfunctions are statistically independent 
and 

(ii) at small energy separations uj ^ 5 ~ {I'oL'^)^^ the result should be of the same order 
as for the one single wavefunction (w = 0). 

This immediately gives rise to the relationship between aq and 7g given in Ea. p7p . 
In particular for q = I one finds: 

71=7 = 1-^- (38) 

The scaling relationship Eg. ([55)) has been first suggested by Chalker [7l| and checked 
numerically in a number of works (77| . Fig|3] gives an evidence that the scaling relation- 
ship Eg. ([57)) between 7, and aq holds true for q = 2 as well. An important point is that 
finite-size corrections are quite large and should be taken into account for an accurate 
determination of the exponents 7^ (an example is presented in Figs l5l6l below). 

So far we have discussed the correlation functions containing |V'i(r)P and thus inde- 
pendent of the random eigenfunction phase. The simplest phase-dependent correlation 
function is involved in the density-density correlation function: 

X V.(r)^*(r')V'*(r)V',(r')) (39) 



According to Ref. TlJ the Fourier-transform of this correlation function is equal to: 



1 D(u},q)q'^ , , 

which is a generalization of the correlation function in the diffusion approximation for 
the case where the diffusion coefficient may depend on q and w. The main assumption 



7l| here is that at the critical point and for q <ti 



-1. 



i.(.,,)=,-F(,L.) = {^4|::u^ £>>!■ (41) 

One can easily see that this assumption is equivalent to the statement that at distances 
i < \r — r'\ < the correlation function K{uj,r — r') has the same r — r'-dependence 
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Figure 3: (Color online) Scaling relationship between 72 and for the 3D Anderson model of the 
orthogonal symmetry class at two different system sizes L=14 and L=24. The theoretical slope computed 
from Eqs.ll3fll,(l26]l gives 72 fa 1.14 ±0.03. 



Eg. (P?]) (but with L — > as explained above) as the correlation function of the local 
DoS Ci(cj,r — r') that we discussed above. This requires non-trivial phase correlations 
in the eigenfunctions which allow to replace the phase-dependent average by its phase- 
independent counterpart: 

{MrMir'W,{r)^Ar')) ^ (|V'.(r)nV,(r')P) (42) 
The replacement Ea. ip^ illustrates the general rule 

Rule (iii): To estimate the averages that cannot be expressed in terms of \ip\'^ (phase- 
dependent averages) one should find by permutation of space and/or energy variables the 
corresponding phase-independent average and apply rules (i) and (ii) in order to estimate 
the latter. 

One should note, however, that this rule does not apply for |r| > L^:. At such distances 
the phase correlations are no longer present, and the phase-dependent average vanishes 
exponentially while the r-dependence of its phase-independent counterpart saturates. 

2.2.2. Scaling estimates for matrix elements: mobility edge. 

Now let us show how to use phenomenology of multifractal wavefunction statistics to 
estimate the matrix elements of local interaction Mijki given by Eq.®. The simplest 
one Mi = Man, so called super- diagonal matrix element, is proportional to the inverse 
participation ratio P2: 

{Mi) ^ ir^'^^'^^'^ L-^\ (43) 

Here and below we define microscopic length-scale i via its relation to the upper energy 
cutoff E(j = The factor 3 in the above equation is of the same origin as in 

the statistics of eigenvectors of real Gaussian random matrices. Applicability of such a 
relation to the case of fractal wavefunctions was demonstrated in Ref. [lH for the case of 
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Anderson transition on a Cayley tree. Qualitative picture, emerging from their analysis, 
is that typical wavefunctions can be represented as products ipj{r) = Xj(r)<i>^(r) where 
Xj(r) describes fast-oscillating functions specific for each eigenstatc j, but insensitive to 
the vicinity of Anderson transition, whereas $^(r) has a meaning of smooth envelope 
function, which is however sensitive to the proximity of energy ^ to the mobility edge E^- 
Combinatorial factor 3 in Ea. (P5|) is due to the fast-fluctuating "Gaussian" component 

The diagonal matrix element Mij = Mijij with i ^ j can be easily found from 
Eq.dSZl); 

{M,,)^V-'M{co)^V-' (^^y , (44) 

where V ~ L"^ and uj — £ i — . The detailed numerical analysis of this matrix element 
has been done in Ref. \7l\ . Here we present in FiglS]data for Ci(uj) = M{ui) correlator 
for 3D Anderson model with Gaussian disorder [W ~ 4) at the mobility edge, computed 
for three different system sizes. The plot of the density of states for the same model is 
presented in Fig. HI The values of 7 extracted for all three sizes are shown in the inset, 
together with extrapolation to L — > 00, which leads to 

7^0.57 ± 0.02 (45) 

Below we will use this value of 7 in our analysis, as well as the value of the prc-factor 

£;o = 2.08 ± 0.25 (46) 

which is also extracted via large-L extrapolation, as shown in FigjBl In addition, wc 
present in Figl5]data for the same correlation function at _B = 8.0 in the localized part 
of the spectrum (to be discussed below in Sec. 12.2.51 and used in Sec. [51). Here we just 
note that the logarithmic slope of the Af(a;) function (rather deep inside localized band) 
is close to its value 7 for the critical eigenstates. 

A straightforward generalization of Eq. pi| is the matrix element 

d^r(|V'.,(r)p...|V'.„(r)p)^V-("-i) (^^j , (47) 

where 7„ is defined in Ea. p7p . and all energy differences — f ; , | are assumed to be 
of the same order lo. 

The off-diagonal matrix element Mijki with all indices different and all energy dif- 
ferences of the same order uj strongly fluctuates and therefore has zero average value. 
However, its average square modulus can be computed as following. Wc start by writing 

{\M,,ki\^) - J d^vd'v' (48) 

(tA* (r)^. (r')V^; (r)V', {v'Wu (rOV-fe {^Wi (r^V-/ (r)) 

Now we apply the rule [Hi) of the previous subsection and define the phase-independent 
correlation function: 

(|V'.(r)nV',(r')nV'fc(r)n^,(r')P> (49) 
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Figure 4: Density of states for the 3D Anderson model with Gaussian disorder with the width W = 4 
for the system size L = 10 (blaek squares) and L = 20 (blue stars). 




Figure 5: (Color online) Correlation function M{u)) for 3DAM with Guassian disorder and lattice sizes 
L = 10, 14, 20 at the mobility edge E = 5.5 (red squares, blue triangles and black diamonds) and at the 
energy E = 8 inside localized band (green dots). Inset shows 7 values for L = 10.12.14.16.20. 
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£„(oo)= 2.08 + 0.25 
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Figure 6: Pre-factor Eq for 3DAM with Guassian disorder width W = A at the mobility edge Ec = 5.5, 
for lattice sizes L = 10, 12, 14, 16, 20 




Figure 7: (Color online) Scaling with the system size of the off-diagonal matrix element. The theoretical 
exponent is —3/2 according to Ea.l l5H . 
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Figure 8: (Color online) Scaling with energy difference \Ei — Ej | ~ \Ei — ~ \Ei — Ei\ ~ ... ~ a; of 
the off-diagonal (|Mijj.jp) (blue dots) and the product of semi-diagonal (MukiMjjfii) (black squares) 
matrix elements. The theoretical exponent found from Eas. ll51l l. ll52ll and Il26ll is equal to —2.13 it 0.1. 



Using the rules (i), (m), we estimate the phase-independent correlator as: 

(50) 



\ / ]^ \ d—a2 / J \ 2{d—d2) 



V4 V|r-r' 

In order to obtain an estimation for {\Mijki\'^) one has to integrate Eg. ([50)) over r and 
r — r'. The first integration is trivial and results in the factor V. To estimate the result of 
the second integration we note that according to Eg . ([551 in a typical sample the exponent 
d — 0^2^ in the power law in Eq. (|50p is smaller than d. This means that the integral over 
r — r' is dominated by large distances \y — y'\ ^ L^^. Finally we obtain: 



V3 



(|M..fcd'>typ - ^ f - ) . (51) 



The same rules applied to the correlation function {MikuMjijk) of semi-diagonal 
matrix elements give exactly the same answer: 

Td / p \ 27 



(M,HiM,,,fe)typ -i^yi;) ■ (52) 

However, the square modulus of the semi-diagonal matrix elements is much larger: 

(|M,fe,z|')typ - Ltr^^-'^^V-^-^-'^iE^/uy. (53) 



The relevant phase-dependent correlation function and its phase- independent counterpart 
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Figure 9: (Color online) System size scaling of semi-diagonal matrix clement for the 3D AM. The 
theoretical exponent found from Ea.l l53ll . ll26t is equal to —1.21 ± 0.02 for the simple averaging and 
— 1.24 ± 0.02 for the typical averaging. 




Figure 10: (Color online) Energy difference scaling of semi-diagonal matrix element for the 3D AM with 
two system sizes L = 14 and L = 24. The theoretical exponent found from Eq.| I53II . II26II is equal to 
-1.57 ±0.03. 
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according to rule (iii) are: 



{\Mr)\'Mr)^Hr)mr')\'r,{r')Mr')) ^ 
mrr\MrT\Mr)mir'r) (54) 
The latter for £ |r — r'| can be found in the following form: 

{\Mr)\'mr')?\Mr)\'\Mr')\') - 

^ (f^) (^) iv^l) ■ ^^^^ 

At |r — r'l > we assume that the r — r' dependence in the first factor in r.h.s. of 
Eg. ([55]) saturates. 

Indeed, according to rule (i), at |r — r'| ~ L S> the following decoupling can be 
done: 



V\mrr\Mr')\'\Mr)\'\MrT) ^ 

V'mrr\Mrr)mr')f\Mr')\') — ] ■ (56) 



^0 



This suggests that V2 ~ '2,{d ~ ■ 

At |r — r'l ^ i according to rule (ii) all wavefunctions are effectively in one space 
point: 

V4(|V,(r)nV;,(r')n^fc(r)n^,(r')n 

« v\m^)\'\My^)?m-r)?) 

t) (^) (") 

Comparing Eqs.(l571),(l5ni) with Eq.(l55l) we find: 

1^1=3(^2-^4), = 2(rf - (i2), v-i^d-d2- (58) 

The |r — r'l dependence in Ea. ((55|) is the power law |r — tc'\~^^^^^^'> with the exponent 
[d — Q2). As in the case of the off-diagonal matrix elements considered above, for a 
typical realization of disorder this exponent is smaller than d. Therefore, the integral 
of the phase-dependent correlation function over (r — r') which determines the matrix 
element (|Mifciip)typ is dominated by large distances |r — r'| ^ (for |r — r'| > L^^ the 
phase-dependent correlation function decays exponentially). Thus plugging |r — r'| 
in Ea. (j55p and multiplying the result by and V (for integration over the remaining 
coordinate r) we arrive at the announced result Eq. (j53p . 

We conclude our analysis by noting that if instead of typical averaging the full ensem- 
ble averaging is performed, the exponent = 3d4 — 2d2 may be (and in the 3D Anderson 
model of orthogonal symmetry class is) negative, due to the contribution of rare untypi- 
cal realizations. In this case the principle contributions to all correlations functions come 
from small distances |r — r'| ^ I. This changes the estimates in Eqs.([5T |) . ([5^ . ([5^ : all 
those equations will acquire an extra factor 

[ujlE^r-l^ (59) 
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which explicitly depends on the fractal dimensionality di and not only on d2. Note, 
however, that in the particular case of the 3D orthogonal symmetry class, the exponent 
|Q!2/rf| ^0.1 is extremely small, so that an extra factor Eg. ((5^ is of order one for most 
of the practical purposes. 

The predictions made on the basis of rules (i)-(iii) of the algebra of multifractal states 
arc checked by numerical diagonalization of the 3d Anderson model of the orthogonal 
symmetry class and summarized in figures FigiTl-Fig llOl One can sec a very satisfactory 
agreement for exponents of various power laws which were found numerically and derived 
theoretically using only one fractal dimension d2. 

Closing this subsection we conclude that the rules (i) — (iii) and the definition of 
multifractal dimensions Ea. (P5| constitute the full set of rules necessary to estimate any 
correlation function of critical wavefunctions. This set of algebraic operations will be the 
main analytical tool to deal with the strongly disordered case considered in this paper. 

2.2.3. Scaling estimates for matrix elements: multifractal insulator. 

Estimates for the matrix elements on the localized side of the Anderson transition can 
be obtained from the corresponding formulae of the preceding subsection provided that 
the localization radius is larger than the characteristic length i (multifractal insulator). 
The only modification to be done in all the local averages (...) is to replace L by Lioc 
and to add a factor (Lioc/i)"'' (where n is the number of different eigenfunctions in 
the average) that accounts for the probability for a point r to be inside the localization 
radius of each of the wavefunctions. To apply these simple rules for the matrix elements of 
different wavefunctions with the energy separations a;, one has also to make sure that all 
the corresponding length scales = [vqU))^^/'^ are smaller than the localization radius 
Lioc- This sets up an important condition which determines the frequency domain of the 
multifractal correlations {multifractal frequency domain): 

Eo > > 6l ^ [lyoLLr'. (60) 

In particular, under this condition we have: 

(A/,) - 3r LfJ% (61) 

which is much larger than the corresponding critical value Ea. (|43p . The combinatorial 
factor 3 in Eq. arises because of the statistics of phase of wave function as explained 
above, see note after Ea. (|15|) . 

However, the estimate for the average diagonal matrix element Mij does not change 
and is still described by Eq.([33]). The reason is the additional (compared to Mi) factor 
{L\oc/L)'^ due to a small probability for two wavefunctions being localized close in space 
which exactly compensates for the replacement of V — > Lf^^ in Eq. (jU) . One can easily 
check that the estimates for the average square of the off-diagonal matrix element Ea. ([5T|) . 
and the average of different semi-diagonal matrix element (|52p also remain unchanged 
under the condition Ea. (j60p . 

The average square of the semi-diagonal matrix element gets enhanced with respect 
to the critical case Eq. ([55|) : 

{\MU')tyv> - ^^"^'"'^^ L-ot Lt (EoM-'. (62) 
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This enhancement factor of {L/L\ocY^ resuhs in a regular scahng of {\Mikii\'^) cx 
with the total volume. 

Note, however, that in all cases the typical value of a matrix element |M'^p| of well 
overlapping states is much greater than the typical average value (|M|)typ, and it can 
be obtained from the corresponding expression for (M)typ (or (M^)typ) by replacing 
the total volume V by the localized volume L^^^. In particular, the typical value of the 
diagonal matrix clement for well overlapping states is: 

^r-Aoc (f^)'- (63) 

This difference between the average value of a matrix element and the typical value for 
well overlapping states is due to the fact that the most of matrix elements in insulator 
are very small due to poor overlap of the corresponding states. 

One can easily check that a,i u = 5l the typical absolute value of matrix elements 
of well overlapping states i, j, k, I does not depend on the number of different states and 
has a order of magnitude of the inverse participation ratio: |Af*^P| ^ £^'^ {£/ Lioc)'^^ ■ 
For uj > Sl the matrix elements M*^^ of well overlapping states get suppressed, and 
the suppression is stronger when the number of different states in the matrix element 
increases. 



2.2.4- Scaling estimates for matrix elements: multifractal metal. 

On the metal side of the Anderson transition, the wave function is globally not a 
fractal (or multifractal), as the moments Pq are proportional to L~'^'^'^~^\ However, 
the correlations of different eigenfunctions (with the energy difference lo) show the same 
power-law w-behavior as the critical eigenfunctions [t^] provided that £ <^ L,^ ^ -Z^corr- 
The physical meaning of the correlation length Lcorr is a typical size of a fractal element 
of which the entire eigenfunction support is built. Thus, locally a wavefunction in a 
multifractal metal is identical to the one in a multifractal insulator inside the localization 
radius. However, the global normalization J^r l''/'(i")P = 1 requires the reduction of 
by a factor of L'^^^.^./V compared to the case of insulator. Thus, we can formulate the 
rule (valid provided that £ <C <^ ^con ) for the estimation of the matrix elements in 
the multifractal metal if their counterparts in the multifractal insulator are known. One 
has (i) to multiply the result for insulator by a factor {L'^^^.^./V)'^ , where q is the total 
number of \ip\'^ in the matrix element in order to take into account the normalization 
and (ii) to eliminate the overlap probability factor {Lf^^.^/V)'^ , where n is the number of 
different wavef unctions, which is no longer needed for the extended metal states. Thus 
the overall factor to add is {L^^^^/V)'^~^ . This immediately leads to 

{M,)^^{L,,,,/£f~''\ (64) 

imkul'hyp - ^ {L.orr/ir-"''^ Lt {E, / LO^ . (65) 

The averages (Mij), {\Mijki\'^), and {MikuMjijk) do not change and have the same order 
of magnitude on the both side of the transition provided that the condition Ea. ((60)) (with 
Leon- replacing Ljoc) is respected. 
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2.2.5. Matrix elements of the off-critical states beyond the multifractal frequency domain. 



As we have seen in the previous subsections, in the multifractal metal and insulator 
characterized by the large correlation/localization length Liocicon- ^ ^, the wavefunc- 
tion correlations are very similar to those of the critical multifractal states at the Ander- 
son transition point. However, this correspondence is only valid if the energy separation 
Lu between the states lies in the multifractal frequency domain, 5^ < uj < Eq, bounded 
by effective level spacing, Sl, and high frequency cutoff, Eq . If one or several energy sep- 
arations are beyond the multifractal frequency domain, the frequency-dependent factors 
in the estimates change. If all energy separations are larger than Eq the frequency- 
dependent factors decrease very fast, so Eq provides the high-frequency cut-off for all 
wave function correlations. 

More delicate is the case where all energy separations are smaller than 5l. The 
behavior of the matrix elements in this region is different in the multifractal metal and 
in the insulator. In the multifractal metal, the w-dependent factors simply saturate fr^ - 
In contrast, in the multifractal insulator they acquire additional powers of ln((5i/a;). 

The logarithmic enhancement factor appearing in the insulator reflects the Mott's 
physics of the resonance mixing of states; this effect is also responsible for the logarithmic 
factors in the expression for the low-frequency conductivity: a '-^ uj In^'+^^i/w). The 
key element of this phenomenon is that localized states with close energies cannot be 
considered independent even when the distance between centers of localization is large 
compared to iioc- To understand the origin of the correlations, we repeat the Dyson 
arguments for the level statistics. Consider two typical localized states with small energy 
separation 6E < Sl at distance R > L\oc from each other and vary the disorder potential. 
The states cease to become orthogonal to each other, their typical overlap decreases as 
t 5l exp(— i?/Lioc) with distance. The energy splitting of the two states becomes 
Lo ~ {(SE)^ + i^)^/^. If i? < L\oc ln((5L/'^) the typical overlap between two states at 
distance R is larger than to which implies that small energy splittings are dominated by 
rare events when t ~ SE ~ w. This implies that in a typical situation the states at 
distance R hybridize forming superpositions i/'±(r) = cosa±V'i(r) + sinQ!±?/'2(i") of the 
parent states ?/'i,2(r) with a± ~ 1. This resonance hybridization makes even remote in 
space parent states mix with each other, provided that the distance between their centers 
of localization is less than the optimal one i? < i?o ~ iioc 3> ^loc- 

The definition of Rq convenient for numerical study is given in terms of the dipole- 
moment matrix element: 



The definition (j66l67D is useful in the range of relatively well- localized states, at w ^ Jl- 
One can easily check that for the extreme strong localization |'0i,2(r)P = ^(r — ^1,2) the 
matrix element Xij = xi — X2 corresponding to ipij = ip± is equal to the x-component of 
the distance between centers of localization of the parent states. The results of numerical 
computation of (|S7)) for the 3D Anderson model are presented in Fig l2.2.5l for Ep = 



Rliuj) = d 



(66) 



where 




(67) 
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Figure 11: Logarithmic dependence of Rq on ui for the 3D Anderson model with W = A and Ep = 8.0 
well in the localized region. At this value of the energy the localization radius Li^c = {'^0^l)~^^^ ~ 3 
(this estimate follows from the data for shown in Fig. |4] and for 5^ shown in Fig. I27II . The linear in 
\nuj fit at small corresponds to R{oj) = AIiiSi^/llJ + B with A = 1.9±0.05, B = — 8.1±0.4. It is valid 
in the region where ij(aj) is larger than 2Lio(,. 

8.0. As expected, i?o('^) is linear in In(cj) at very low values oicu/SL- However, there 
is a broad transient regime with essentially non-linear in In a; behavior for moderately 
small Lo. It is important to note that the values of R{lo) in this regime are smaller or 
comparable than 2Lioc (estimated from the inverse participation ratio, with the use of 
data for the DoS value and typical level spacing Sl ) which is the minimum distance 
between centers of localization where the Mott's physics of resonance mixing strictly 
applies. 

It was shown in Ref.jT^I that the matrix element My is proportional to i?g^^(a;). 
This result can be also justified by the perturbation theory arguments similar to those 
used in the derivation Rq ^ L\oc ln((5L/aj). Consider the parent states that are strongly 
localized at sites m,n and have a distance R = |Rm — R-nl between centers of their 
localization. We now vary the realization of disorder and treat the matrix element of the 
corresponding change of the disorder potential Hmn as a perturbation. The amplitude 
of the cigcnfunction |'i/;m(Rn)P at a (remote from its center of localization Rm) site R„ 
becomes of the order of 

The matrix element is given by 
M^JO^ E El'/'™WI'l^"WI'^2|V'™(Rn)poc ('^') r dRR^-'e-''/'^^-^. 

parent states 7' ^ ^ 



Here instead of summing over the parent states we integrate over the distance between 
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the centers of localization R taking into account the statistical repulsion of centers of 
localization at i? < i?o ~ -^loc ln((5L/a;). The final expression in terms of the integral 
over R is similar to the corresponding expression for the Mott's frequency-dependent 
conductivity but differs from it by an extra R^ because of the square of the dipole 
moment matrix clement and an extra in front of the integral. The estimate of the 
integral in the limit Rq ^ iioc finally gives Mij oc Rq^^ and cr(uj) cx lo"^ ^o^^- 

Summarizing the results of this analysis we conclude that the Ea. (|44|) for Mij in the 
multifractal insulator becomes 



-5 J ^»Sl 



{M.,)=V-' (^) { ^'-J' ^ (68) 



A simple analytic expression that smoothly interpolates between these two asymptotes 
in Ea. (p5)) can be written for the function M{uj) = V{Mij) with w = — ^ji 

M{u:) = , (69) 

ln'(^+c) + l 



where the constant c > 1. One should also take into account small variations of the 
effective fractal dimension d!^ as Ep moves away from the mobility edge or the disorder 
parameter W moves away from the critical value [77| . This results in the dependence of 
the effective exponent: 



= l-d2/d + a{5L/Eo), a>0. (70) 

The energy scale 6l can be expressed through the inverse participation ratio {Mi) in the 
form which is convenient for numerical simulations: 

Finally, the upper fractality scale Eq can be found from the condition V {Mij) w 1 at 

w = So. „„„ 

The enhancement, Eas. (|68l69p . of the overlap of two localized wavefunctions results 
from the fact that the eigenfunctions which are anomalously close to each other in the 
energy space (closer than 5l for the localized wavefunctions) are automatically well 
overlapping, as they arc just the symmetric and anti-symmetric combinations of the one 
and the same pair of localized wavefunctions with the optimal distance between their 
centers of localization being R^^ ^ L\oc ln((Si/a;). The most important consequence of 
this phenomenon is that in 3D case the correlation function M{uj) continues to grow at 
low frequencies as log^ ^ . which makes it possible to establish the superconductive order 
in some part of the domain Tc <C (5l , as will be discussed below in section [21 



3. Insulating state. 



In this section we discuss the physical properties of the insulating state in the vicin- 
ity of superconductor-insulator transition and show that it is characterized by the large 
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single particle gap which is responsible for the activation temperature dependence of 
conductivity, (7(T) oc exp(— Tj/T), observed in many works (46. 30l 21 at low tempera- 



tures on the insulating side of the transition. We start by assuming that Fermi energy 
Ep is deep inside the region of localized states, so that the interaction of electrons from 
different localized orbitals is weak and leads merely to a small perturbation. Then the 
Hamiltonian © can be further simplified to 

where the scaling estimate for the typical value of matrix elements Mj = J dr^j{r) is 
given by Ea. ((6T|) above. We will refer to the last term in (|72|) as to the local pairing 
coupling] formally it looks like the "negative-U" local attraction considered in Q. 

The operator product in the last term in Eq. (j72p is identical to the occupation number 
product Uj-fUji which is equal to 1 if both available electrons states are populated, and 
to otherwise. Thus the only role of the interaction term in (j72p is to shift down energies 
of all double-occupied orbitals. Note that one does not encounter such a term in usual 
theory of disordered superconductors, since Mi vanishes in the thermodynamic limit for 
delocalized electronic states (compare Eqs. (|5T|) . ((M| ). 

Let us order all eigenstates i/jj according to the increase of eigenvalues . Then the 
last filled eigenstate ipm of the Fermi-sea (at T = 0) for even total number of electrons 
is defined by inequality 

2Cm - gM„^ < < 2^™+! - gM„,+i (73) 

(we count all single-particle energies from the Fermi-energy). In a macroscopic system, 
the energy interval in ([75]) vanishes as the inverse volume, 1/V. Within the "even sub- 
space" of the whole Hilbcrt space (i.e. each orbital is either empty of double-occupied), 
the local pairing can be fully accounted for by the redefinition = — f Af^. 

However, single-occupied orbitals are not involved in this interaction. 

The increase of thermodynamic potential ft due to addition of odd electron to the 
ground-state is 

(5r2oc = = £,m+l — 6ri+l + 6ri+l = (74) 

|m™+i + o(v-1) 

Using the Eq. (pT|) we estimate typical value of Sftoc- 

sniyj' = ^g£-%Uoc/e)-''\ (75) 

where i/oc is the localization length for the states with E = Ep. 

Consider now two single-electron excitations on top of a fully paired Fermi-sea defined 
by Ea.(|73p. which are produced by transferring of one electron from the m-th state to 
the m + 1-th one. The energy of the two single-particle excitations (which results from 
the depairing) is 

2A(r^ = U+1 ~U+ gMm = |(M™ + AWi) + O(V-i) (76) 
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Figure 12: (Color online) Distribution of the inverse participation ratios P2 for the 3d Anderson model 
at the Fermi energy Ep = 9.0 on the insulating side of SIT (the mobility edge Ec = 5.5) for different 
system sizes. 

thus the typical value of the pairing energy Ap is also given by Eg. ([75)) : 

AP = Igi-HLUer'^ = yi^o {^^y^ (77) 

In the right-hand side of (j77p we reintroduced the dimensionless coupling A = gvQ and 
the parameter Eq = l/(i/o^'^) which determines the high-energy cutoff for fractal corre- 
lations and we also employed Ea. p^ . The average single-particle density of states z^(e) 
is determined by the probability distribution V{M) of inverse participation ratios Mj. 
Namely, the probability to find an excitation with the energy e coincides (in the limit 
V — > 00 ) with the probability to find a value of inverse participation ratio M < 
Therefore 

iy{e) = va V{M)dM (78) 

Jo 

We have generated the distribution function V{M) numerically, using the three dimen- 
sional Anderson model with the Gaussian distribution of local energies, Ea. ((23l) with 
W = 4. The mobility edge in such a model is located at \E\ = 5.5. 

Numerical data for the distribution function 'P{M) are shown in FiglT^lfor several 
sizes L and the Fermi energy = 9.0 in the localized part of the spectrum. These data 
demonstrate a sharp drop of V{M) at the values of M much smaller than the typical 
value M'yP, as well as a considerable size-dependence of the slope. FigfT2l shows that the 
low M tail can be well approximated by the exponential dependence 

V{M) (X e-'=(*-f'"/^^) for M < A/'^p (79) 

where the coefficient c depends on the energy Ep. We used the finite-size scaling together 
with an extrapolation to large L limit in order to get the values of c appropriate for a 
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Figure 13: (Color online) The exponent c of the small-IPR tail of the distribution function P(M / M^^^) 
at the energies Ep = 8.0, and 9.0 for different sizes L, together with extrapolated value Coo for each 
energy. 



macroscopic system. The results obtained for energies Ep = 8.0 and 9.0 are shown in 
FiglTSl The data shown in Fig 113] demonstrate that the extrapolated (at L — >■ oo) value 
of the coefficient c saturates in the interval 

csa 0.73 ± 0.05 (80) 

for the Fermi energies deep enough in the localized band. This is in agreement with the 
one-parameter character of the distribution V{M) = V {M / M^^^) in the L ^ oo limit. 

In Fig. [T3]we also present the linear scale data for the distribution function P{y) of 
local gaps Ap,(at the Fcrmi-cnergy Ep = 8.0), normalized to the typical value Ap^. Note 
that the data both for large and for small y can be fitted quite accurately by an analytical 
expression shown on the plot which contains exponential factors and l/y^ dependence 
relevant for the intermediate 1 < y < 6. This modification leads to a somewhat lager 
coeflicient c in the exponential dependence e~'^l'^ . 

We use the numerical data for V{M) to obtain, according to Ea. ((78)) . the average 
density of single-electron states. The results are shown (for the two Fermi-level positions 
Ep = 8.0 and Ep = 9.0) in FiglTS] Both plots for DoS i^(e) nearly coincide while 
expressed as functions of the reduced variables e/Ap in the most interesting region 
X < 1.5. At the lowest energies DoS decays exponentially, i^<<(e) ~ VQe~'^^''/'^ . The 
similar results for the low energy DoS were obtained in [5[ (see Fig. 16b of Ref. Q). 

Practically, the shape shown in Fig ll5l implies the existence of a nearly hard gap 

^^^c,^p^^\E^(^l—El.^ ' ci«0.2-0.3 (81) 

The average density of states could be measured directly by the tunneling conductance 
via a large-area tunnel junction. The problem with such measurement in an insulator 
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P(y) = {^//)exp(-c/y-b,y-fc/) 
A= 1.80 ± 0.03, c= 1.08 ±0.01, 
fa = 0.15 ± 0.02, 6, = 0.030 ± 0.003 
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Figure 14: (Color online) Distribution function for the normalized local gaps y - 
energy Ep = 8.0 on the insulating side of SIT (the mobility edge Ec = 5.5). 
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Figure 15: (Color online) Gap in the one-particle density of states due to local paring of electrons on 
the localized states at the energies Ep = 8.0 (black stars) and Ep = 9.0 (blue circles). Scaling with 
X = Mi/Mtyp for relatively small x is demonstrated. 
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is that electrons should be evacuated somewhere after tunneling, otherwise the strong 
Coulomb blockade would make the measurement impossible. A possible method to avoid 
the Coulomb blockade might be to study tunneling conductance through a relatively thin 
film with both its surfaces weakly coupled to metal contacts. 

The local tunneling conductance measured by STM is expected to show a threshold 
behavior with A^'' corresponding to the state m localized near the observation point 

in real space and with close to the Fermi energy. The local gap A^-* fluctuates 
from point to point and is distributed according to Fig llSI We are not aware of such 
measurements in the insulating state, the data on the superconducting side of SIT are 
given in Ref. [42[, for temperatures both above and below Tc- We present detailed results 
for the tunneling conductance as a function of temperature in section 16.21 

Above the transition temperature the data [42| show large suppression of the density 
of states without any coherence peak that appear at the gap edges below Tc. The 
absence of coherence peaks at the gap threshold shows the qualitative difference between 
the local pairing gap due to simple binding of two localized electrons and the BCS gap 
appearing due to the many-body correlations within the energy range Abcs around 
the Fermi-level that one expects in the small grain of superconductor. In the case of 
fractal localized states all filled (double-occupied) levels are shifted down in energy, thus 
the total number of states near the Fermi-energy (in a stripe of several Ap width) is 
not conserved and there is no reason for appearance of a peak above the gap. In that 
sense, the local pairing (due to attraction!) plays the role similar to that of the Coulomb 
repulsion in suppressing the tunneling conductance. 

The above results were obtained neglecting all matrix elements of interaction except 
from super-diagonal ones, Mj. We now discuss the validity of this approximation. Indeed, 
in the problem of ultra-small grains treated by Matveev and Larkin [s^l, the diagonal 
matrix elements Mij were taken into account while calculating the parity gap (fT4|) . via 
the renormalization (jl3p of the local pairing energy. This renormalization was necessary 
(even for A <C 1) due to the ultraviolet (UV) divergency in the Cooper loop diagrams. 
Pairing on the fractal eigenstates is of different nature: the frequency-dependence (|^^ of 
matrix elements Afy eliminates the UV divergence, thus the virtual i — >■ j transitions can 
be neglected and the high-energy cutoff w/j is not necessary as long as 6l ^ T^- Here 
T° is the superconducting transition temperature for Ep at the mobility edge Ec given 
(within modified mean-field approximation developed below in section 2]) by Eq. (|164p . 
Using of Eqs. (|71l81ll64p one can eliminate the model parameters Eq and A and express 
the local gap Ai via the observable quantities and 6l'- 

_ 3ci 6l ^ Sl , . 

where for the 3D Anderson transition 7 = 0.57 and C = C7(0.57) ~ 3.1. The relation ([5^ 
shows that in the insulating region 5l 3> due to the nonzero fractal exponent 7 the 
local gap value Ai grows with 5l much weaker than 5l itself. In the limit 7 — > Eq. (|5^ 
transforms into the analog of the Matveev-Larkin relation ([Ti)) , after a sub- leading term in 
the denominator of (|82p is taken into account via a substitution 7 {5l /T^V ^ i^L /T^^V ^ 
1. 

Now we turn to the discussion of the intrinsic low-temperature conductivity of the 
insulator with localized pairs. The binding of electrons into local pairs diminishes the 
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single-particle DoS and thus suppresses the variable-range hopping conductivity. A clas- 
sical example of such an effect is the " Coulomb gap" due to Efros and Shklovsky: the 
soft gap i'{E) ~ in the average DoS leads to a transformation of the Mott law 
a{T) cx e-^'^"/'^^^'^ into the Efros-Shklovsky law a{T) cx e-('^«s/^)''^ In our case the 
low-energy states are exponentially rare, and their account leads to the logarithmic in 
temperature corrections to the activation energy determined by the hard gap ((8T|l : 



a{T) (X exp 

where 



Tln{To/T) 



(83) 



To = |Ap « 0.25Ap (84) 

Note that the nearly activated behavior given by Eas. (|55|) and ([M]) can hardly be distin- 
guished from the purely activate one with the activation gap Ai given by ([ST|) or (|5^ in 
the limited temperature range available in most experiments. 

We thus associate the spectral gap Ai with the measured [i^,!!^! activation energy Tj. 
The external parameter {Ec — Ep) representing the disorder strength in Eq. ((8T|) can be 
replaced with an experimentally more accessible parameter (1 — a /uc) oc [Ec — Ep)/ Eq. 
Here a is the high temperature conductivity and CTc is the value of the conductivity where 
the parity gap Ap first develops. Wc obtain 

Ti = A{1- aja^y^^ , A k. 0.5XEq (85) 

where A is conductivity-independent. This equation predicts a moderate increase of T/ 
with disorder strength in agreement with the experimental data (46| . see Fig. 1161 The 




0~(kFl)(au 



Figure 16: Experimental values of the gap from Ref . 1461 . T/ (boxes) and a fit to the equation II85II with 
u^v = 1, d2 = 1.3. 

only fitting parameter used in Fig. [TBI was the value of the constant A = 0.5Ai?o ~ lOK. 
Assumin g th e applicability of the BCS theory for (less disordered) superconductive InOx 
samples [46[, and using [83[ the estimates wd ~ 500 A' for Debye frequency, we find 
A « 0.2 and Eq ~ WOK. 

Applicability of the scaling formulas to the above analysis is not obvious. Indeed, 
the value of 1 — (t/ctc ~ 2/3 for the most insulating sample on the plot of FiglT6l We 
have demonstrated, however, that the correlation function M{uj) (see Fig which is 
closely related to the inverse participation ratio, is approximately described by the critical 
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scaling even deeply inside the insulating phase. Then assuming that Ea. ((T6)) also holds 
approximately in a relatively wide interval of 1 — cr/cc (although, perhaps, with some 
"average" exponent i^av instead of ly) we obtain an extrapolation of Eg. ([55)1 for an entire 
interval of ct/ctc relevant for the experiment. Note that the classical experiments on Si:P 
system have shown that unless the special care was taken to shrink the interval of l — a/uc 
to be much smaller than 1, the observed scaling of localization radius was characterized 
by the exponent i^av < 1-5, sec Fig. [2] That is why it is not surprising that reasonable 
fit to the experimental data shown in Fig. [12] corresponds to the choice u^^ = 1 . 

There are a few reasons why the estimated value of w lOOA' turned out to be low 
in comparison with the Fermi energy Ep ^ O.SeV for amorphous InOa; samples with the 
electron density in the range of 10^^cm~^. This value implies that ratio Eo/Ep ~ 0.03 is 
one order of magnitude lower than our estimate (|46)) obtained for the 3D Anderson model 
with Gaussian disorder (which, however, should not be expected to be quantitatively 
correct for a-InOa;). First, this estimate of Eq uses BCS formula to related the value 
of the critical temperature in the less disordered samples Tc ~ 3K [i^ to the value 
of the coupling constant in all materials of these series. The BCS formula is however 
qualitatively wrong for the samples close to mobility edge, which are in the regime 
discussed in the section |4l In this regime the critical temperature is given by (|164p which 
leads to a much lower value of the interaction constant A « 0.02. The estimate for the 
upper energy cutoff Eq becomes lOOOA' which is in agreement with (^5)) . This estimate 
neglects the effect of thermal fluctuations discussed in section [4.2.51 the effect of these 
fluctuations is to reduce Tc so the actual value of the interaction constant corresponding to 
the sample close to the mobility edge with Tc « 3K might be slightly larger A « 0.03—0.05 
corresponding to Eq « 400 — 6QQK. These values of Eq/Ep ^ 0.1 — 0.3 are roughly what 
one expects from the analysis of the three dimensional Anderson model. Finally, wc note 
that appearance of low energy scales in the insulating samples of InOa; was conjectured 
in the early paper (s^ for completely different reasons. All these arguments demonstrate 
that the values of A and Eq that wc obtain from the fit of the experimental data are 
roughly what one should expect in these samples. 

The reasonable fit to the data was made possible by a small value of the exponent in 
Eg. (1851) w hich is substantially less than lyd 3 due to eigenfunction fractality. Thus, the 
data [46[ provide the indirect evidence for the eigenfunctions fractality in the insulating 
samples of InOa,. 



4. Cooper instability near the mobility edge: the formalism. 

In this section wc develop approximation schemes to treat the Cooper instability and 
superconducting order formation in the regime when the Anderson theorem 3^, 3^ is 



not valid due to very strong disorder. Namely, in section 14.11 we give two versions of 
the modified mean-field approximations, MFA, for determination of the superconductive 
transition temperature T^, in section [4.21 we derive an analog of the Ginzburg-Landau 
functional which is necessary to estimate the role of fluctuations beyond the modified 
MFA, and in section 14.41 we give an alternative method of Tc determination: the virial 
expansion. The nature of approximations involved in the modified MFA and the virial 
expansion methods are very different, thus reasonable agreement between the results 
obtained by these methods indicates the validity of both of them. We show that in a wide 
region near the mobility edge the dependence of the superconducting critical temperature 
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on the interaction constant can be found analytically with the accuracy up to a pre-factor 
of the order of one despite of the presence of strong thermal and mesoscopic fluctuations. 
In a more disordered sample the critical temperature can be determined using the semi- 
analytical approach of virial expansion applied to the pseudo-spin Hamiltonian. The 
results of this section allow to suggest the pseudo-spin Hamiltonian that provides a unified 
description of superconductivity in homogeneously disordered (non- granular) systems 
including the BCS regime and the region close to the superconductor-insulator transition. 

4.I. Modified mean-field approximation. 

The goal of this subsection is to develop two versions of a modified mean-field ap- 
proximations that can be used to determine the critical temperature of a fractal super- 



conductor. We start from the standard Abrikosov-Gor'kov- Anderson [35|, |36| mean-field 
equation for Tc of a disordered superconductor: 

A(r) = J d'^v' K(r, r') A(r'). (86) 

The kernel K{r, r') is defined by 

Kiry) = I ^Ty,, ^,(r)^,(r)V',(r')V,:(r'), (87) 

where '(/'j (r) arc exact single-electron wavefunctions and 

tanh(e;/2r) + tanh(C,-/2r) _ 

fl%3 = T-TT ' = ^^^1 

K'l + Kj 

Within the standard mean-field approximation, one neglects the spatial variations of 
A(r). Then integrating over r in Ea. (|5B|) and using the orthogonality and normalization 
of different single-particle wavefunctions one can eliminate all ^Jiiv) out of Ea. (|5B)) . Note 
that it is only possible if the wave functions are real, V'* (r) = V'i (i") ; the time-reversal 
symmetry is preserved. This leads immediately to the well-known equation for in 
terms of DoS function i'(^): 

1 = 1 E^^^^/ ^M^MO (89) 

i 

where j;, = 77,1 ~ ^'^^ tanh(^i/2r). 

However, the approximation of a constant A cannot be used under strong disorder 
conditions near the mobility edge. Physically it is due to strong mesoscopic fiuctua- 
tions of the local DoS function i^(^,r), see Below wc propose a modification of the 
MFA scheme which makes it possible to account for the major part of mesoscopic DoS 
fiuctuations. 

To construct the modified MFA, we note that the appearance of the solution to the 
Eq.dHH) is equivalent to the divergence of the series Tr(l — K)^^ = X^i^Lo^r^'"- By 
d'Alambert criterion the latter is equivalent to the condition 

lim — — - 1. (90) 

n— >oo irA 
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Explicitly, the trace Trif" can be written as 

TrX" = (.9/2)" (91) 

where Mij^i is defined in Eq.(l8]). 

Now we make a crucial approximation: we neglect all the off -diagonal matrix elements 
in Eq.(|S]) with more than two different indices. The diagonal matrix elements My = 
Miijj = Mijij = Mijji with only two different indices will be retained. By so doing, 
we neglect thermal fluctuations of the superconductive order parameter; we also neglect 
a part of mesoscopic fluctuations, whereas the most important type of the them (i.e. 
fluctuations of the local density of states) will be taken into account. Thus our approach 
has a meaning of a modified mean-field approximation. The neglected off-diagonal terms 
determine the strength of fluctuations on top of the modified MFA solution. We derive 
the corresponding Ginzburg-Landau functional and estimate the role of such fluctuations 
in subsection 14.21 below. 

There arc two ways to impose the "diagonal constraints" the for matrix elements 
Mijki- (i) to set ii — ji thus imposing n constrains in Ea.()92p. or (ii) to set ji = ji' or 
ji = ill with I 7^ V . In this case one imposes 2n — 2 constrains, so that only 2 summations 
will remain. As any summation gives a macroscopically large number of terms, one can 
neglect all terms corresponding to (ii) in the thermodynamic limit N ^ oo. Thus, upon 
neglect of off-diagonal matrix elements one finds eventually 

Tri^"«(5/2)" ^ 77,^77,,. ..77,„ M,„,,M,,,,...Af,„_,,„, (92) 

ii...in 

where Mik are defined in Ea. lfTU)) . It is easy to show now that the condition given by 
Ea.(|90)). with Trif" defined by Ea. (|92)) . is equivalent to the solvability condition for the 
equation 

A,; = I ^ Afc r^fc Mfe„ (93) 
fe 

which is the basis of the modified MFA we will be using below. 

The set of new order parameters A^, entering Eq. (j93|) . represent the superconducting 
"ordering field" acting onto a pair of electrons occupying the i-th orbital ^pi{r). The idea 
of the approximation is that, instead of using a constant (in real space) order parameter 
A(r) = A, we assume a smooth dependence of A^ on the single-electron energies ^i. 
Under such an assumption, together with the replacement of Mik matrix elements by 
their averages (according to Eg. ip^ ). the Ea. (P51) can be transformed into the integral 
equation 

A(0 = ^ I dC^(C)M(e ~ C)A(C) (94) 

Eg. dM]) is a natural generalization of the BCS mean field equation which follows from 
it at M{uj) — 1. However, allowing for the energy dependence of M{uj) may lead to 
drastic consequences. Indeed, a simple scaling analysis of this equation with "critical" 
M{uj) = {Eq/uj)"^ given by Ea. (|44)) leads to the power-law dependence of the critical 
temperature on the attraction interaction constant A: 

T, (X EoX^^^ , (95) 
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where 7 is defined in Eas. (|38l45|) . 

We encounter an unexpectedly strong increase of T,. in the small A limit with respect 
to the usual result, Tc oc exp(— 1/A) for a conventional BCS superconductor (with the 
same A). The price to pay for this increase is a very strong inhomogcneity in the real 
space of the local pairing amplitude: 

A(r) = |5]Afe,yfcV^(r). (96) 

k 

The corresponding analysis will be presented in subsection 15.31 

To demonstrate the validity of the solution one needs to plug it into Eq.(|55)) 

and to use "diagonal approximation" 

y"d''r^;V,(r)V'^(r) = A/,i<5y 

together with Ea. ((93)) . 

The transformation from Ea. (|M|) to Ea. (|M)) is not exact: we replaced the fluctuating 
matrix elements Mik by their averages according to Eq. pi]) . Thus Ea. (|M)) contains an 
additional (with respect to Ea. (|93|) mean- field- type approximation. One can eliminate 
this additional approximation. For this purpose, let us define the new matrix Q and 
vectors 4>i: 

Qik = ^y/ViVkMik (t)i = /^iy/rfi (97) 

Then the solvability condition for Ea. (|M)) transforms into the condition that the largest 
eigenvalue fcmax of the symmetric matrix Q becomes equal to unity. This condition is 
equivalent to the instability onset with respect to formation of a superconducting order 
parameter A^, as will be seen below in subsection 14. 2. II 

In the thermodynamic limit the condition fcmax = 1 implies global superconductivity 
only if the corresponding eigenvector of the matrix Q is extended. Otherwise, fcmax = 1 
means formation of local " islands" of new phase with uncorrelated local order parameters 
in these islands, it does not immediately result in any global order parameter (for the 
discussion of relevant examples see [s^, H^). Eigenfunctions of the Q matrix may 
be localized, if the effective "coordination number" (the bandwidth) Z of this matrix 
stays finite in thermodynamic limit. General estimate for the coordination number is 
Z ^ v^TcLf^^ because the relevant energy window populated by "active" single-particle 
states ( taking part in formation of the many-body superconductive state) is of the order 
of Tc, and each eigenstate | i > is coupled by the the matrix elements Mij to all neighbors 
\ i > in the localization volume Lf^^. This argument shows that as long as the single 
electron eigenstates are delocalized, Ep < Ec and within diagonal approximation^ the 
effective coordination number Z = 00 and no localization of the Q eigenstates is possible. 
However, this conclusion becomes invalid when the off-diagonal matrix elements Mijki 
are taken into account: the lowest eigenstates of the full kernel K{r, r') get localized. The 
most important physical consequence of this localization is that, as the temperature is 
decreased, the superconductivity first appears in small, well separated regions, similarly 
to the situation realized in superconductors with inliomogeneous Tc (87| . In the regime 
of delocalized single electron states a further decrease of temperature results in a more 
homogeneous superconductivity. As the Fermi energy is increased past the mobility edge, 
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Figure 17: (Color online) Density of states p{k) of the matrix Q, with peak near k = 1, for L = 14 (black 
squares), L = 19 (blue dots) and L = 22 (red diamonds). 

single electron states get localized and the issue of the Q-eigenfunction localization and 
resulting inhomogeneity become relevant even within the " diagonal approximation" , we 
discuss this regime in section [51 

Here we present, as an example, the numerically obtained spectrum of the Q matrix 
for a finite system at the mobility edge. Fig. [T7] shows the averaged (over 2000 real- 
izations) density of states for the Q matrix generated for the 3D Anderson model with 
Gaussian on-site disorder of the strength W = A, at the energy E = 5.5 corresponding to 
the mobility edge. Three system sizes, L = 14, 19, 22 were analyzed. In a finite system 
the peaks of the finite width are seen in the density of states p{k) near fc = 1. One sees 
that with the system size increase, the peaks become more narrow, so it is natural to 
assume that they evolve into a J-function peak in the L —¥ oo limit. The corresponding 
temperature (chosen so that to pin the DoS peak position to A: = 1) is thus associated 
with Tc. 

In subsection 15.21 we compare results for the critical temperature Tc obtained by 
the three methods: modified analytical MFA equation numerical generation of the 
spectrum of the temperature-dependent matrix Q, with the temperature T = Tc adjusted 
so to locate the peak in p(k) at fc = 1, and virial expansion method, to be described in 
subsection 14.41 . 

Our modified MFA method is similar in spirit to the usual BCS mean-field approx- 
imation: it neglects both thermal and mesoscopic fluctuations of the pairing field A(r) 
with respect to its " background configuration" . The essence of the modification we used 
is that our background configuration is not uniform in real space. Instead it is given by 
the function A(r) defined in Ea.ljM]). In order to estimate the strength of fiuctuations 
we neglected, we will develop, in the next subsection l4.2l a Ginzburg-Landau description 
for the order parameter field configurations A(r) which are close to the background field 
A(r). 
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4.2. Ginzhurg - Landau functional. 

Here we derive the free energy functional _F[\E'(r)] defined in terms of a smooth enve- 
lope function 

*(r) = A(r)/A(r) (98) 

and valid near the transition temperature Tc for ^'(r) field which smoothly varies in 
space. 

We start by writing the action in terms of the Grassmann fields (p„ 

5 = ^ ip*^{-iu>n + Hi)ipn + ^ tp%^{iuJn + -ffl)<^-„ (99) 
n n 

"l,n2,n3 

where Hi is the single-particle part of the Hamiltonian, given e.g. by Ea. ((20)) or Ea. ((2T|) : 
g is the interaction constant and n and ni^2,3 are the Matsubara frequency summation 
indices. We decouple the interaction term via the Hubbard-Stratonovich transformation 
involving the auxiliary field A(r; r2„). Below we focus on the static n = component of 
this field. In doing so we neglect quantum fluctuations of the order parameter but take 
into account the thermal and mesoscopic fluctuations. 

We expand the Grassmann fields over the cigcnstates '0i(r) of the single-electron 
problem 

3 

and represent the thermodynamic potential S'[A, $] as 

^ r n ij ^ ^ ^ 

where 

= ^(^) V'.(r)^,(r) = Y DkRk^J. (102) 



and 



Here we defined 

Rk^J = j dri^fc(r)^,(r)'/', (r) Dk = J drA(r)7Afe(r) (103) 

Note that the way we decoupled the interaction term is specific to the superconduc- 
tive correlations. A generic Hubbard-Stratonovich field should contain also a compo- 
nent V{r;iln) coupled to the combination (x!"')*x!"'' which corresponds to interaction 
in the particle-hole channels. Wc neglect such interactions in this subsection and will 
study their effect in subsection 16 . 1 1 in connection with the 5,f5| terms in the pseudo-spin 
Hamiltonian. 
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The matrix variables A^j , which appeared naturally in the second term of the action 
are not mutually independent. Indeed, the number of independent components of the 
(discretized) field A(r) scales with system volume as cx V, whereas the number of matrix 
elements scales as V^. The matrix elements Ay arc mutually constrained due to the 
ortho- normality conditions / dripi{r)ipj{r) =Sij. 

Now one can complete the Hubbard- Stratonovich transformation by performing the 
Gaussian integration over to obtain the Tr In of the corresponding matrix in Eq. (|10ip . 
In addition, one can represent the first term in Eq. ljlOip in terms of the coefficients Di 
of expansion of A(r) over the full set of single particle wave functions ipi{r). Expanding 
the Tr In over Di up to fourth order, we obtain free energy functional in the form 

= -J2 D*ICuD, + \Y1 D^Dj ^^^^^ ^fe^*' (104) 



0^ ' ■' 4 

ij ijkl 



where 



In Eas. pIH)) .([TU5 |) we denote 

Jijkl = ^ ^ ^iViVlRjuiV2^kU2V3^ll^3Vi CyiV2,V3Vi (106) 



where the function ('(■••) is defined by 



CviV2,V3.'^i — Cv2Vl:'^il'3 — Cv3l'i,l'iV2 — (10"'') 

77 rrz r VlviV2 + ^1^31/4 ~ 1^3 f 2 ~ Vi^ii^i] ■ 

Below we will use Eas. (|104p - (|107p to derive the effective Ginzburg-Landau functional 
^[^^(r)] in the form 

FGi[*(r)] = J^oT^Jdr (^a{r)'f^{r) + ^^\r) + C\W'i>{r)\^^ (108) 

The factor i^oT^ in (|108p is introduced so that to keep ^E" and a(T) dimensionless. The 
physical properties of the superconductor described by functional (|108p arc controlled by 
four parameters: a = (a(r)), b, C, and by the strength W of fiuctuations 6a{r) — a{r) — a 
defined as 



W = J dr {5a{r)Sa{r')) = (109) 

= V-^ J dvidv2dv3dvi{5K(v4,V2)5K{v3,V4)). 

where (5/\ (r, r') = A'(r, r') — (A'(r, r')). We are able to describe the spatial fluctuation by 
one parameter, W , because (5a(r) correlations are short-ranged compared to the typical 
scale of ^'(r) variations. 
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4-2.1. Transition temperature: coefficient cl(T). 

We begin by evaluating the coefficient a that vanishes at the transition point, a(T) = 
a{T — Tc). For this calculation it is sufficient to use a constant 'l'(r) = vp. Now let us 
transform the first term of Ea. (|104[) . using Eqs. (jl05l96ll02ll03l) . Wc assume here that 
Aj obey the matrix equation (j93p . i.e. the local order parameter coincides with A(r) 
defined in Ea. ([M)) . The result reads 



g ^ g2 ^ 



(110) 



I J ijk 

where we used the " fusion rule" following from the completeness of the set of tpj (r) : 

RkiiRpqi ~ Mkipq. (Ill) 



and an equivalent form of Ea. (|93p written in terms of Di and Ajk variables: 

A = f Xr;fciAfc,i?,fez. (112) 



kl 



To determine the coefficient a we neglect the off-diagonal terms Muki with fc 7^ / in 
Ea. (|110p and thus reduce it to 



(113) 



where Q and (j>i are defined in Eg. ([77)) . In the continuum limit one writes <^,; — ^ ^, 
Mij — > pA^(^ — C)i Si ~^ J d£, and thus the bilinear form pi3p transforms into 



(0) _ i^oV 



2 / / d^dCLi^xmm) 



where 



with the kernel 



A 



(114) 

(115) 
(116) 



The function 4>{£,) is, by construction, an eigenfunction of the integral equation with the 
kernel Q{£,,C)'i its eigenvalue k(T) approaches 1 at T = Tc. At this stage we need to 
specify the normalization condition for the function 



(117) 



The form of the normalization condition Eq. ()117p was chosen to enable a smooth crossover 
to the conventional (non-fractal) BCS case 7^0, when (f>j=o = Tyjrji^. To get the 
equations valid in the crossover regime we need to introduce the high-energy cutoff Q.d so 
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that all integrals over rf^ go over the range \^\ < fi^. Then a straightforward modification 
of the normalization condition 



„ In ^ 



-0.0 + \ 



(that reduces to (|117p when 7 ^ 1/ln^ = A) allows one to recover the BCS limit at 
at a fixed and 7 — > 0. We are mainly interested in the case of intermediate to strong 
fractality, 7 ~ 0.6, and weak interaction A <C 1, so we shall use the simplified condition 
(|117p appropriate in this limit. 

Making use of Eq. (|117p wc reduce the expression (|114[) to 

^ = VQVTH^-i~^[k{T) - fc2(T)] (119) 

Since 1 — k{T) ^ 1 at T « Tc, it is sufficient to evaluate the derivative dk/dT\k=i- We 
find it by rewriting the kernel (|116p in a dimensionless form; 



2\x~yp\j xy 

where At = \{Eq/2T)'^. Clearly, dlnfc(r)/(iln(r) = -7, and, finally, the coefficient 
a(T) is expressed only in terms of the ratio of T /Tc'. 

a{T) - (120) 



Note that the order parameter ^ introduced in Eg. ([98]) is dimensionless, whereas A has 
a dimension of energy. 



4-2.2. Quartic term: coefficient b 

Now we turn to the computation of the coefficient b in the expansion (|108p . Our 
starting points are Eas. (|104p and (|106ll07p . As in the previous part we may use here 
'l'(r) = const, thus A(r) is proportional to A(r) defined in Ea. ([M|) . Therefore using 
Eas. ((96|) . (jl03p we may substitute Di = ^ J^a ^aVaRiaa h^ the fourth order in Dj term 
of Eq. (|104p reducing this term to the following form: 

^4 = \[^yj2Y.Y.^-^'>AcAdVaVbVcVd (121) 

ijkl abed i^i ■■1^4 

Now we do summations over i,j,k,l using the fusion rule (jllip . and then neglect the 
off-diagonal matrix elements, i.e. we set everywhere M^-^^^aa <^i/ii/2-^fia- Summations 
over a, b, c, d can be done now with the use of modified MFA equations ((93|) . Finally, we 
proceed from the last remaining summation over vi to integration, and find 

52 



where 



n3 1- /- 1 / tanh(a;) 1 



<^{x) = 16T^ lim = -J 



with the function Cijfc; defined in Ea. (|107p . In tlie limit 7 — > we have A(^) = const and 
the integral in (|122p reduces to 

- / dx<^{x) = ^ - ifeBCS , 

327 ^ ' 167r2 2 

leading to the standard BCS result for the h coefficient in Eq. (|108p . In general case we 
find 

b=\J dx x^ coth^(a;) $^{x) (123) 

where 0(x) = T^^/^(f>(2Tx) and (/)(^) is the fc = 1 eigenfunction of the kernel (|116p 
subject to the normalization condition (|117p . 

4.2.3. Gradient term: coefficient C. 

We now turn to calculation of the gradient term of the functional (jlOSp . Note that 
previously employed "diagonal approximation" for the matrix elements Mijki is not suf- 
ficient for that purpose. In order to find the coefficient C we need to take into account 
the off-diagonal matrix elements Mijki which contain pairs of different levels with nearby 
energies: j^j ~ Cil ~ ICfc ~ 61 ^ with LL!{q — >■ 0) — )■ 0. The reason for that is similar 

to the one which applies to the standard BCS Hamiltonian with the reduced interaction 
term g p' c-f.p c^.-p c]^ p' c| _p, . One has to allow for the non-zero momentum of a 
pair and thus go beyond the "diagonal approximation" in p and p' in order to be able 
to compute the phase rigidity C which is related with the supercurrent. 

We first illustrate this statement using the standard theory of disordered supercon- 



ductors as an example, see e.g. [88|. In that case, quadratic term of the free energy 



expansion (jl04p can be represented in the coordinate space 

i(5(r-r') - A'(r,r') 
.5 



F[A(r)] = j dvdv' 



A(r)A(r') (124) 



where 

Kiry) = .^Y ^^^U^^^^^ ,125) 

Here a;„ = T:T{2n + 1) and (|125p is the equivalent form of Eqs. (|87l88p . As the stan- 
dard procedure goes, one averages the kernel (jl25p over disorder in the semiclassical 
approximation and obtains the Fourier-transformed kernel 

^('^)=2.-oTE^^,^ (126) 

where D is the diffusion coefficient. Finally one obtains for the coefficient C in Ea. (|108p 
the well known result 

TT D 



dq^ 
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q^O 



In order to generalize this derivation to the case of strong disorder in is convenient to 
introduce two particle spectral function 



X {Mr)4'Jir)r^{r')rA^')) S{Ef~C^). (128) 
The average kernel in Ea. (jl25p can be expressed through this two particle function: 

K{q)^vaTy I — : , , . r 129 

(4 -iw„)(f + 

In the diffusive limit the spectral function is 

1 Dq^ 
'^"^'^^^nn^ + iDq^r- 

and its substitution into Ea. (|129p gives back Eq. (|126p . 

This derivation of the standard result demonstrates that the small g-dependcncc of 
the averaged kernel K{q) comes from matrix elements of the operator e*''^''"'" ' between 
eigenstates ■!/'i(r), '0i(r) with nearby energies Cjl ^ Dq^ . Notice that diagonal approx- 
imation used in previous sections corresponds to z = j in Ea. (|125p . In this approximation 
one neglects the off-diagonal terms Mijki completely and gets g{q, Vt) = S{il)f{q), where 

f{q)^ I dvdv' {i,l{v)^l{v')e^'^^^-^">) 



For delocalized wavefunctions which occupy all the available volume one can write 

mr)^K^')) « (^.^(r)) {^j{v')) = 1 + 0(l/.9), 

where g 3> 1 is the dimensionless conductance Then f{q = 0) = 1 + 0{1/ g) and f{q ^ 
0) = 0(1/5). The jump in f{q) at g = implies that in the diagonal approximation the 
coefficient C is infinite in a metal which is obviously wrong. However, the value of K at 
q — Q can be determined correctly in the diagonal approximation because at q = the 
spectral function g{q, Vt) — > (5(51) both in the diagonal and in the diffusive approximation 

Ea.dnnn. 

The correct derivation of the coefficient C in case of a strong disorder has to take 
into account a non- uniformity of the "background configuration" of the superconducting 
order parameter given by Eq. ([Ml) . Wc follow the same logics as in the derivation of the 
a and h coefficients in previous sections. Note that the g-dependence comes from the 
second term of the matrix (|105p which leads to the second sum in Eq. ()110p . We focus 



on this term and repeat the steps of the derivation which led to Ea. ljllOp using now a 
weakly modulated 'l'(r) = ^^e*'^'' with small q. The second term of ()110|) becomes (in 
the continuous form, after replacing — > Vv{)d£^): 

j dciF^^\q) = - J dql^g]^^ J d^dujduJiduj2A{C++u;i)Ai(^+u;2) (131) 
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where ^± = C i '^/2; the function A(^) = 4>{0 / VviO obeys Ea. (|94| and is normalized 
according to Ea. ()117p . and the function M4{q;Lu,LUi,uj2) is defined by 

\ijkl I 

where = ii- 

h,kM) = J drrfr'e^'i^'-'"-') Mr)i'A^)Mr>Ar')^lirHfir') (132) 

The function A/4(g; wi, aj2) replaces function g(q,uj) employed in the standard 
derivation of the jV^E'p term for diffusive superconductors. The additional terms ip'^{r) 
and ipf(j') (and the summation over the corresponding states k and I) arise due to the 
essentially r — dependent " background configuration" which r-dependence is expressed 
through ■)Jj'^{r) by Eq.dM]). 

The main source of fiuctuations of M4 are the spatial fiuctuation of the eigenfunctions 
whereas the spectral fluctuations are relatively weak. Neglecting spectral fluctuations, 
one can switch from summation over indices k, I to integration over the energy vari- 
ables ^i, ^j, ^fe, ^; with the constant DoS and consider the average /(q, {w}) = (lijuil)) 
as a function of all energy differences uj, wi, uj2- Then the expression for the coefficient 
C in the Ginzburg-Landau functional Eq. (|108p takes the form 

C = -^^^^0 j d^dujdujiduj2^{^ + (^/2 + uji)A{^-uj/2 + uj2) (133) 
XTj{^ + uj/2, C - uj/2) Tj{^ + uj/2 + uji) Tj{^ - uj/2 + L02) J({c^}), 
where J({c.}) = -V^ ^^^^ . 

^ q=0 

The next step is to estimate /({w}) using the rules formulated in section 12.21 Wc 
begin by considering the critical wavefunction statistics. Applying the rule (iii) wc find 
the phase-independent counterpart to the combination of eigenfunctions in Ea. p32p . It 
happens to coincide with the one for the square of the off-diagonal matrix element Mijki 
and is given by Eas.(|l ^ . ([5n)) . Then expanding Ea. p32p up to and integrating over 
(r — r') up to |r — r'l ~ L^, we obtain an expression similar to Ea. ((5T|) but containing an 
extra factor due to the expansion of e''^^''^'" up to q^: 

ji{u})^v'Lim,,,w„^Lt+^ ^^^^^ 

In Ea. (|134p we assumed all energy differences to be of the same order wi ~ a;2 ~ The 
estimate Eq. (|134p holds true on the insulator side of the Anderson transition as long as 
w ^ <^L (see the discussion in section [2.2.41 and Ea. (|5(I)) ). 

Now we have to estimate the result of integration over the energy variables in Ea. (jl33| ). 
Using the asymptotic behavior of A(cl>) cx oj"^ which follows from ((94|) . (see section |4T|) 
and ri{oj) oc w"^ and Eq. (|134p power counting shows that for Tc > Sl the integral is 
dominated hy u ^ uji ^ uj2 ^ Tc and can be estimated by 
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Remarkably, the enhancement factor depending on the fractal exponent 7 is canceled by 
the coupling constant A due to E g. (| 95 p . and the result is not sensitive to fractality and 
is essentially the same as in Refs.[4|, |3[. This cancellation occurs due to the presence of 
the additional (with respect to the standard expression Eq. I128P factors -0^ (r) and ij^f (r) 
in the integral Iijki{q) defined in Ea. (jl32p : without these factors, the final result for C 
would contain extra small factor ~ (Tc/Eq)''. 

Note that the coefficient C is dominated by the off-diagonal matrix elements only 
in metal or in very weak insulator. As one moves towards a strong insulator where 
Tc < 5l, the main contribution to C becomes the one that originates from the diagonal 
approximation and can be roughly estimated as 

C ^ lL, {5l > T,) (136) 

in the diagonal approximation. The estimate (jl36p is based on the simplest picture that 
wavcfunctions localized at distances larger than localization length Lioc do not overlap. 
In fact, as we discuss in section [6.11 this approximation is a bit too crude as it misses 
an important logarithmic factor which increases C values in the range ^ 5l, see 
Ea.(fT751). 

4-2. 4- Mesoscopic fluctuations: coefficient W . 

Disorder always leads to spatial fluctuations of parameters which enter the Ginzburg- 
Landau functional; the major effect is due to fluctuations of a{T, r). Universal mesoscopic 
fluctuations (which provide a lower bound for the strength of this effect) were studied in 
Ref. 0] for usual disordered superconductors and more recently in [3l[ for 2D films with 
the strong Finkelstein effect. Here we follow the same general approach but implement 
it in the Fock space instead of the coordinate space. 

We start from the Eq. piOp for the quadratic part of the free energy. Previously we 
neglected off-diagonal matrix elements Mijkk entering the second term of (jllOp : now 
our goal is to estimate the strength of mesoscopic fiuctuations produced by these matrix 
elements. Thus we represent the second term as 

F22 = E ((^')^^ + '^''^J- (137) 

where Q is defined in Ea.(|Wl) and 

= ^ E ^'^i^' ^'^'^^rj m {mm ) ^ (i38) 

k^l 

The matrix contains corrections from the (previously neglected) off-diagonal matrix 
elements. Its average Tij contributes to the shift of the critical temperature which can 
be estimated using scaling arguments developed in section [2.2.21 The result is that the 
relative shift due to the off-diagonal matrix elements is {5osTc) /T^ ^ 1. 

Here we are interested in the strength of mesoscopic fiuctuations W defined in 
Ea. (|109p and thus calculate the free energy cumulant: 

{v^T^fW - -l-((i^22F22» = (139) 
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Using Ea. (jl38[) and switching to a continuum representation, we obtain the irreducible 
correlation function ((Fr)): 



((^11, 12-^^3,^4)) — 

(t'oV)" y'ViiVi2Vi3'nU J ■■■ J d^md^nd^kdClVninVkl j ... j d'^Vi(f v2d'^V3d'^Vi 

where the double brackets ((...)) denote the cumulant average defined by (we use below 
indices ii... instead of V'lV-- etc.): {{iii2izi4mnkl)) — {i 112131 A'mnkl) — {iii2mn) {i^iikl) 
(1) 

and (r) = '0™(r)?/'„(r). 

We now estimate the coefficient W for the critical states near the Anderson transition 
point. First of all wc note that the decoupled average {iii2mn) {i^i^kl) depends only on 
the two of the three independent differences in coordinates Rs^s' ~ \i^s ~ ^s'l = 
1, ...4), while the cumulant average depends on all the three of them and vanishes when 
Rs,s' > (as before wc assume that all the energy differences are of the same order w). 
In the region where all i?,,..,' < the decoupled average is smaller than the cumulant 
one. Assuming all the differences of coordinates are of the same order Rs^s' ^ R and 
applying the rules (i) — (iii) of section [^.2. II we obtain aX R < L^: 

^l{v,)i:l{v2)^l{v,)^l{v,)^(^^^^^^^ 
^l{v,)i^l{v2)^l{v,)i^l{vA)^^^^^^^^ 

«(ri)V^f,(r2)^^3(r3)^fJr4)^^(ri)^^(r2)V^^(r3)^f(r4))^ 



where a2 = Sd^ — 2d2, = 7ds ~ 6^4 . Indeed, at i? ~ eigenfunctions in dif- 
ferent space points are statistically independent and the corresponding averages can 
be decoupled. The resuh is V"^(L„/^)''('^-''2) in agreement with Eq.([32l). At E ~ ^ 
all eigenfunctions can be considered in one space point and the averaging then gives 
V~^(Ltj/^)^^''~''**^ The last estimate coincides with Eq. (P5|) at small energy separations 
Li^ L where one can consider all eigenfunctions to be identical, and it corresponds 
to all eigenfunctions independently averaged at large energy separations when ^ i. 
The decoupled average estimated in the same way using Eas. (|49)) . (j50|) is of the order 
of V"** {L^ / Rf'^-'^"''' (L^/i)^'^'^-''^') and thus is smaller at i? < than the cumulant 
average by the factor {L^/ RY^°'^ . 

Now we estimate the result of the four spatial integrations: over difference of coor- 
dinates i?i2 = |ri — r2|, i?i3 = |ri — rsj and i?23 = jra — r2| and one free integration 
which results in the factor V. At this point it is important that the dependence of the 
cumulant average on i?i2, i?i3, i?23 is symmetric and such that the power of R is typically 
smaller than Sd, as a^^P > (see Ea. (|55|) ). This means that the main contribution to 
the integrals over i?i2, i?i3, -R23 comes from the region R12 ~ i?i3 '-^ R23 ^ L^. This is a 
crucial circumstance that eliminates the dependence on higher fractal dimensions. 
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Thus we arrive at 



The remaining energy integration is estimated in the same way as in the previous sub- 
section. As before, the dominant contribution comes from the energies within a strip of 
width Tc near the Fermi energy. The final result reads 

((^^22^^22)) ^ A-^ z^oVr^ ^ ~ y^VTl (140) 

Again, as in Ea. (|135p . the enhancement factor which depends on the fractal exponent 7 
cancels out exactly by the coupling constant A, and the final result for the coefficient W 
is 

W^-k^, {5l^T,). (141) 

This result was obtained using the definition (|109|) which makes sense if the spatial scale 
La of a(r) fluctuations is small compared to the scale of 'l'(r) variation. We expect that 
the same estimate p4ip is valid if both length-scales are of the same order. 

In the limit of strong insulator, one can repeat the above analysis to arrive at 

W ^ LL, {5l > Tc). (142) 

This result is obtained in the diagonal approximation. 

4-2.5. Ginzburg parameters for thermal and mesoscopic fluctuations 

Now we use the results given by Eas. (|120ll23ll35ll4ip to estimate the relative width 
of the fluctuation region near the thermal transition into a fractal superconductor state. 
First we estimate the Ginzburg parameter Gi which determines the reduced temperature 
range |1 — T/Tc\ < Gi where thermal fluctuations are strong in a 3D system [89|: 

Gi^^^^^^l, iS,«T.) (143) 

The relative width of "smearing" of the superconductive transition due to positional 
disorder is given by the parameter Gid defined as follows (see e.g. Ref.fSTt): 

Gid - - 1, {Sl < Te) (144) 

The estimates (jl43p and (|144p demonstrate that the modified mean-field approxima- 
tion developed in this section can be used (with relative accuracy of the order of unity) 
in order to estimate Tc of a fractal superconductor. Based on this result we conclude 
that in the region of extended and weakly localized single-particle states with Sl ^ Tc 
the thermal fluctuations of the order parameter phase and the mesoscopic fluctuations 
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of the local Tc which were not taken into account in the modified MFA can at most re- 
duce the global Tc by a factor of the order of one compared to the modified MFA result, 
Eq. (IMI) but can hardly lead to a modification of the functional dependence of Tc on the 
interaction constant A. This is the most important conclusion of this subsection. 

Another conclusion concerns the role of the off-diagonal matrix elements in the region 
of extended and weakly localized single-particle states. These matrix elements are com- 
pletely neglected in modified MFA which nevertheless gives accurate results. They are, 
however, necessary for the correct account of the thermal fluctuations and the local 
mesoscopic fluctuations. The off-diagonal matrix elements also determine the electromag- 
netic response and are necessary for calculation of the critical current, as both properties 
are related with the gradient term in the Ginzburg-Landau functional Ea. (|108p . 

In the region Sl ^ Tc the use of the simplified estimate C ^ Lf^^ would lead to the 
conclusion that the parameters Gi and Gid are large: 

where Z^ff ~ Tc/5l is the effective coordination number to be discussed later on in 
relation with the pseudo-spin Hamiltonian. In fact, more accurate estimate for ZcS 
given by Eq. (|179p in section [01 below shows that both Gi and Gid remain of the order 
of one in the broad range of large 5l/Tc ratios. 

4-.3. Pseudo-spin Hamiltonian. 

As we will see below, there is a sufficiently wide range of parameters on the insulating 
side of the Anderson transition where the superconducting transition temperature Tc is 
of the order of its value ri*^^ right at the Anderson transition point while the paring 
gap Ap introduced in section [3] is much larger than Tc- This means that practically 
the entire region where Tc gradually decreases with increasing disorder falls into this 
pseudo-gap regime. The modified MFA does not work in this regime because Ginzburg 
parameters Gi and Gid are larger than 1. However, the problem can be significantly 
simplified by making use of the large value of the gap A p between states with even and 
odd number of particles occupying any orbital V'i(r). Namely, as the creation of the 
odd state with one particle on an orbital takes a large energy 2Ap to break the pair, 
the transitions between the even (having two or no particles on any orbital) and odd 
states described by the off-diagonal matrix elements can be neglected. Technically, it is 
equivalent to neglecting all the off-diagonal matrix elements and considering the sectors 
of Hilbert space with even and odd states as completely decoupled. Restricting ourselves 
to the low-energy even sector one can rewrite [36| the Hamiltonian Eq.Q in the form of 
a pseudo-spin Hamiltonian with the spin operators 5*^^, Sf acting in the Fock space of 
orbitals \i) = ipii^c) and the diagonal matrix elements Mij playing a role of the coupling 
matrix. The most general form of this Hamiltonian is: 

Hps = ^0(2^1 + 1)- + 
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where the set of operators 



St=^Ai 57 = (5+)t ^^=2(^E4.^--lj (146) 

is equivalent to the set of spin-^ operators S — ^a. Here ^ is random energy distributed 
with density in some interval around 0. 

The Hamiltonian Ea. (|145| ) is the basis of theory of the pseudogap superconductivity 
we will develop in section IB] However, it is valid in any case where the off-diagonal 
matrix elements may be neglected for this or another reason. It was originally suggested 
by Anderson [s^ for a BCS superconductor where the off-diagonal matrix elements are 
small in the parameter Tc/Ep. It is parametrically justified in the pseudo-gap region 
where the energy denominator associated with the even-odd transitions is large since 
Ap ^ Tc- However, the estimate of the Ginzburg parameters Gi ^ Gid ^ 1 (see 
discussion in section [6.11 and (jl79p ) shows that it is also useful for a semi-quantitative 
(up to a factor of order one) determination of the transition temperature of the fractal 
superconductor near the Anderson localization transition. We will use this Hamiltonian 
in order to obtain the phase diagram as a function of (T, Ep) in the entire region spanning 
BCS and the pseudo-gap regime. 

4-.4- Virial expansion method. 

In this subsection we develop a new approach based on the virial expansion method 
applied to the pseudo-spin Hamiltonian Eq. (jl45p an use it to determine the supercon- 
ductive transition temperature and thus the full phase diagram of the disordered super- 
conductor. The approximations implied by this method are completely different from 
the ones of modified MFA. It will turn out to be very useful in order to find the pre- 
cise limit of applicability of modified MFA and determine Tc in the region of localized 
single-particle states where the parameters Gi and Gid can be large. 

Developing our scheme, we follow the approach of Larkin and Khmelnitsky [oot who 
first used the virial expansion method to determine the temperature of magnetic phase 
transition in metallic alloys (see also ) . The idea of this method is to express the free 
energy as a series, where each term contains an exact contribution from a fixed number 
of local variables (e.g. spins for the problem of magnetic impurities): 

Y,F., + Y,{F,,-E,-F,) (147) 

i i>j 

Ejk- F,k+E, + Ej+Ek) + ... 

Here Fi = — TlnTre^^*/-'" is the free energy of a single i-th spin in a field, Fij = 
— nnTre^(^'+^j+^'3''/^ is the exact free energy of two interacting spins, etc. The 
terms in the brackets in the second sum and higher order sums cancel each other for 
large space separations (e.g. jr^ — rj| for the second order term) so that corresponding 
variables become essentially independent. 

When the system is approaching a phase transition, all terms of the virial expansion 
become relevant and an actual calculation of critical singularities becomes impossible. 
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n=l 

i>j>k 



However, the virial expansion method can be used in order to find an approximate location 
of the transition point. Indeed, consider the virial expansion for some susceptibihty 



ciT) = J2xn{T) (148) 



ra=l 

which must diverge at T = Tc. The value of Tc can be determined from the condition 
that Tc corresponds to the limit of convergence of the series Eq. p48D : 

lim ^^^^^ = 1. (149) 

Xn[J-c) 

In practice, an exact analytical calculation of Xn with large n is very cumbersome, such 
calculations being usually limited by the first few terms, n — 1,2, 3.... In the following 
we will use instead of Eq. p49p the approximate truncated criterion 

X2{Tc) = X3{Tc). (150) 

This step constitutes the key approximation of the virial method of calculation of Tc. It 
consists in an extrapolation (in general, uncontrolled) into the thermodynamic limit of 
properties found with exact treatment of few-spin systems. 

We now apply this idea to the specific problem of the calculation of superconducting 
transition temperature corresponding to the Hamiltonian Ea. (jl45p . The relevant suscep- 
tibility is defined with respect to the "ordering" field A which enters the Hamiltonian 
via the source term 

Fa = -^(A5+ + A*5-) (151) 
j 

Transition into a superconducting state is signaled by the divergence of the Cooper 
susceptibility 

^(^) = "aA9A^ = ^^"(^^ 

n—1 

Below we will describe calculation of the lowest-order virial expansion terms Xn (T) with 
n = 1, 2, 3 defined in accordance with Eq. (|147p as: 

XI = Ex!'^ (153) 

i 

X2 XI (X nm Xn Xm / 

n>m 

X3 ^ ^ (XnmZ ^7il ^ml Xnm ^~ Xn Xm X/ ) 

n>m>l 

where Xi^\xij\xi^h-- Cooper susceptibilities of the system of 1, 2, 3 ... spins. 

These susceptibihties can be expressed through the eigenvalues Aa(A) of the A^-spin 
Hamiltonian: 

.W = -^Ee-i"VT^^^ (154) 



X' r, 
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where Zq = ^ and 7q is the sensitivity of the a-th eigenvalue of the N-spin 

Hamiltonian to the A perturbation: 



A„(A) = AW + |Ap7„ + o(|An. 



(155) 



In order to compute the A^-spin susceptibihty xi™'; at one has to represent the iV-spin 
Hamihonian as a 2^ x 2^ matrix. The single-spin susceptibility can be found easily: 



ttT 



(156) 



where Ei, the the upper energy cutoff and C = 0.577... is the Euler constant. 
The two-spin Hamiltonian is given by the 4 by 4 matrix 



V 



(2) 



2^12 



A' 
A' 




V, 



A 



A 



+ 2 '-^12 ('^12)* 
2 "^12 



2 ' 

-'12 + yyi 



\ 




A 
A 



(157) 



where ^± = ^1 ± ^2, J-. 



and the perturbation term y cx A. 



Using the Hamiltonian (|157p we calculate A^'^'' and 7q, and then use Eq. (|154p to obtain 

second virial term X2{T). We give here its simplest form corresponding to jj'2 = 0, the 
full answer being given in the Appendix A: 



X2 



J- 



2E- 



2E+E- 



tanh 



tanh 

E- 



T 



tanh 



^tanh 



E_ 

k 

T 



+ 



2E4 



tanh 



(158) 



1 



tanh ^ 



where 



0)^ 



(159) 



are the exact energies of the two-spin problem. 

The three-spin susceptibility X3 requires the solution of the cubic equations which 
is too long to be written here. The corresponding derivation is given in Appendix 
[Appendix A| in a form suitable for numerical calculations. 

For the standard BCS problem with My = y — ?> 0, the renormalized energies E± 
coincide with the bare ones and , and the first term in Ea. (|158p is the leading one. In 
this case the summations over the two energy variables become independent and one gets 

X2iT) = g'^uo In^ (^^^) • Applying the simplest truncated criterion x'^' (TJ = x^^^ (^c), 



(ie^ /'K)Ehe In the same way one can show that 
Xn+i(rc) gives the same correct result for Tc of the 



one finds the correct result: Tc : 
any truncated criterion Xn(T'c) 
BCS problem. 

This is no longer true for the case where Mij has essential energy dependence, as 
it has in the fractal case where Mij is given by Eq. pi)) . In this situation the double 
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sum in xi and the corresponding multiple sums in the higher-order susceptibilities Xn 
do not factorize. Moreover, all Xn with n > 2 are dominated by the low-energy region 
^ ^ T and do not contain any logarithmic divergence, in contrast to xi- This is why the 
simplest truncated criterion Xii^) = X2{T) docs not make much sense: the corresponding 
solution depends on the upper energy cutoff Eb, while for any other truncated criterion 
Xn{T) ~ Xn+i{T) with n > 2 the solution is independent of Ei,. It is for this reason 
that we used Ea. (|150p to find an approximate transition temperature. Note, however, 
that the temperature ri_2 found numerically for the 3D Anderson model with critical 
disorder from the simplest truncated criterion xi{Ti-2) = X2(T'i-2) gives values of 
which are close to T2-3 found from Ea. (|150p for a meaningful values of the upper energy 
cutoff Eb- 

So far we completely neglected the presence of the odd sector of the Hilbert space. 
This is justified in the region of the pseudo-gap superconductivity where the paring gap 
Ap 3> Tc- However, even when being decoupled from the even sector, the odd sector 
cannot be ignored completely in a general case of Ap ^ Tc- The reason for that can 
be easily seen from Ea. (|154p which is valid in a general case too. Decoupling of the 
two sectors manifests itself in vanishing of all "partial" susceptibilities 7q, corresponding 
to the states a of the odd sector. Thus only the eigenvalues Aq°'' corresponding to the 
even sector appear in the sum in Ea. (|154p . This is not true for the partition function 
Zq, where all states contribute. For Ap ^ T one can neglect e~^°'^'^^^ compared to 
g-Aovo„/T g^j^j^ arrive at the result given in Appendix [Appendix A However, in the BCS 



and the fractal region where Ap < Tc the susceptibility will be substantially decreased 
by "statistical dilution" of the even states by the odd ones. 

From the standard theory of the BCS superconductivity it is well known that in all 
expressions the temperature dependence appears as the factor tanh (^^^ ■ In contrast 

to it, the Eas. (|156l I158P contain factor tanh ( |^ ) . This gives a simple rule to adjust the 



virial expansion scheme to the case of vanishing Ap: 

E- E- 
tanh — tanh — (160) 

Thus the temperature T2-3 found from the formalism described in Appcndix |Appendix A| 
should be simply divided by 2 to find the corresponding approximation for the transition 
temperature in the region of extended single-particle states. In Appendix [Appendix B| 
we present the derivation of the virial susceptibilities for the general case of an arbitrary 
paring gap Ap (still, neglecting off-diagonal matrix elements, which mix the "even" and 
"odd" subspaces of the Hilbert space). In particular, we show there that as long as 
the renormalization of energy levels (i.e. the difference between E+ and in Ea. (|159p ) 
can be neglected, one can find (see Ea. (|B.2p ') the whole family of distribution functions 
F(e, Ap) that interpolate between the two limiting cases, Ap ^ Tc and Ap <C Tc. 

In the region of extended and weakly localized single-particle states the renormal- 
ization of energy levels is indeed negligible which justifies the approach made in the 
derivation of (jB.2p in this regime. This renormalization is not small in the region of 
a strong pseudo-gap but in this regime one can project the pseudo-spin Hamiltonian 
Ea. (|145p into a truncated Hilbert space that consists of the even subspace only. Thus, 



63 



both these limiting cases should be well described by a simple replacement rule: 



tanh(|;) ^F(e,Ap) 



sinh (l^) 



(161) 



cosh(f) +6-^^-/"^' 



where Ap = (Ap). Application of this replacement rule to all the formulae obtained in 
the pseudo-spin model is expected to take correctly into account the effect of " dilution" by 
the odd sector. With this modification the pseudo-spin Hamiltonian Ea. (jl45p becomes a 
unified tool to describe (semi-quantitatively or better) the entire region spanning weakly 
disordered BCS superconductor to the Anderson insulator. 

5. Superconducting state very close to mobility edge. 

In this section we derive the main physical properties of the three dimensional dis- 
ordered superconductor with the Cooper attraction characterized by the single-electron 
states that are critical or very weakly localized. As explained above, a very important 
property of these wave functions is their fractality which is responsible for the unusual 
properties of these superconductors. In the regime considered here the number of paired 
electrons in localization volume is still large, which allows us to neglect the level spacing 
Sl compared to other relevant energy scales: 



The most important quantity that characterizes a superconductor is its transition 
temperature and the single particle gap. To compute them we analyze the modified 
MFA developed in subsection 14. 1 1 and determine the critical temperature Tc and the gap 
function A(^). Then in subsection l5.2l we compare Tc obtained by the modified MFA with 



analysis of the discrete modified MF equation with the fluctuating kernel Q, and the virial 
expansion method. We show that the results for the transition temperature Tc obtained 
by different methods are in a good mutual agreement but exhibit a strongly enhancement 
of Tc compared to the expectations based on Anderson theorem. Finally, in subsections 
I5.m5.5l we use Ginzburg-Landau functional derived in subsection l4.2l in order to estimate 
the fluctuation corrections to these results, to calculate the distribution function of the 
local order parameter near Tc, the local density of states at low temperatures and the 
supcrfluid density. 

5.1. Pairing in the modified mean-field approximation. 

The linearized modified mean-field equation for fractal superconductor was derived 
in subsection 14.11 see Ea. ((94|) . It can be generalized for an arbitrary T < Tc where a 
non-zero gap function A(^) develops: 




the results of two other methods based on the " diagonal" approximation" : the spectral 



- C)A(C) 



tanh 



(162) 



2T 



A straightforward calculation shows that near the transition temperature Eq. (|162p leads 
to the same result as the Ginzburg-Landau functional with coefficients found in section 
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Figure 18: Numerical prcfactor C'{-y) appearing in the equation for the transition temperature 1 1164 I I 
obtained by numerical solution of linearized integral equation II163I . 



14.21 above. The properties of the superconducting phase in the vicinity of the transition 
temperature can be found from the solution of the hnearized integral equation, equivalent 
to dMl): 

A«)^^/^c "«-f^'^' .a„l4 (163, 

Eas. (jl62ll63p with the power-law kernel M{uj) given by Ea. ipH) can be solved nu- 
merically. Its solution depends on the fractal exponent 7 which controls the power of the 
kernel. In particular, the critical temperature is found to be 

T,"(A,7) -SoAi/^C(7) (164) 

where the function C{'y) is plotted in FigHS] As was already mentioned above, the 
power-law dependence on the interaction constant A, instead of the exponential depen- 
dence e~ A in the standard BCS theory, implies a dramatic enhancement of Tc by disorder 
if interaction constant A is small. This result could be important for observation of super- 
conductivity in a system of cold fermions with weak attraction trapped in a disordered 
lattice. The solution for A(^) at T = Tc is shown in FiglTSk for 7 = 0.57 corresponding 
to the 3D Anderson model. At large ^ 3> Tc this function decays as A(^) oc as 
demonstrated in Fig|19b. 

The maximum value of the function A(^) is attained at ^ = for all temperatures. 
In the r — > limit we find 

Ao(0) - EoX'^''D{-f) (165) 

where the function -D(7) is plotted in Figl^Dk. The ratio 2D{'y)/C{j), which characterizes 
the ratio of low-temperature spectral gap to the transition temperature, is shown as 
function of 7 in Fig|20b. Coincidentally, this ratio for the 3D Anderson model with 
the fractal exponent 7 = 0.57 turns out to be rather close to the BCS value 3.5. This 
implies that unfortunately it is difficult to distinguish the fractal 3D superconductors 
from conventional BCS ones by measuring 2A/Tc value. 

Evident strange feature of D{j) and 2D{'-f)/C{'j) behavior is that they do not seem to 
approach the BCS limit at small 7. The reason is the same as was discussed in section 
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Figure 20: (a) Prefactor -D(7) in the equation II165II for the maximum value A(0) of the T = obtained 
from the numerical solution of the non-linear gap equation 116211 . (b) Ratio 2A(0)/Tc that follows from 
II165II with D{')) shown here and 11164 l l with C(7) shown in Fig llSI 



14.2.11 around Eq. ([118^ : in order to get the correct crossover to the non- fractal hmit 
7 — > , one needs to introduce the upper energy cutoff ^d- Indeed, the BCS hmit is 
reached at 7 ^ ln~^ <^ 1 which is never satisfied if the upper energy cutoff is infinite 
as in the calculations here. 

5.2. Comparison of Tc values obtained in three different approximations. 

General modified MFA equation (jl64p . with the parameter C{'y) « 3.1 evaluated for 
7 = 0.57, and the estimate (|46| for Eq, gives 

T^" = (6.5 ± 0.8)Ai '^^ (166) 

As we explained above, there are two different types of corrections to this formula: the 
corrections caused by the off-diagonal matrix elements in the original Hamiltonian, and 
the corrections due to approximations made within the "diagonal approximation" , in 
particular, the continuum approximation that neglects the fluctuation of spectrum and 
the matrix elements which are taken to be a power-law function of the energy difference. 
An estimate of the contribution from off-diagonal matrix elements is given in section l4.2.5l 
In order to clarify the role of the static fluctuations of the spectrum, we determine the 
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Figure 21: Transition temperatures determined via p{k) spectrum, for a number of coupling constants 
and system sizes between L = 10, 12, 14 and 19. 

transition temperature for the same 3D Anderson model using two other methods. We 
first determine it from the condition Amax = 1 (see section l¥?T|) where Amax is the largest 
eigenvalue of the matrix kernel Q defined by Ea. (|97p that contains matrix elements 
Mij corresponding to the specific realization of disorder in finite size samples. This 
method adapts modified MFA to the case of strong spatial fluctuations. Second, wc 
apply the virial expansion method as described in section [4.41 Because here we consider 
the regime Ap ^ Tc, wc use the replacement rule Ea. (|160D to take account of the dilution 
of pseudospin states. In both methods we observe strong finite-size effects as shown in 



The extrapolation to infinite system sizes gives the values of Tc as function of the 
dimensionless coupling constant A shown in Figl^Sl As one can see, all these methods 
give the results close to each other and to the result Ea. (jl66[) of the modified MFA 
at 7 = 0.57. The conclusion that we draw from this apparent coincidence is that all 
approximations involved in these methods are reasonable. Specifically, it seems that one 
can neglect the fluctuations of the single-particle DoS and the matrix elements My at 
the mobility edge; one can also use the small-n truncated criterion for Tc in the virial 
expansion method. The main correction to the results of the analytical modifled MFA 
given by Eas. p62p comes from the off-diagonal matrix elements. However, even this 
correction is able to reduce the transition temperature at most by a factor of the order 
of unity leaving all the results of modified MFA semi- quantitatively correct. 

5.3. Pairing amplitude in the real space 

The real-space pairing amplitude corresponding to the solution of Ea.([Ml) can 
be determined from Eq. (P5)) for T ^ Tc. To demonstrate that A(r) corresponds to 
strongly spatially inhomogeneous solution, we calculate the averaged square of the pairing 



FigEnn 
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Figure 22: Transition temperatures determined via tlie virial expansion, for a number of coupling con- 
stants and system sizes between L = 10 and L = 22. 




Figure 23: (Color online) The results for transition temperature as function of dimensionless coupling: 
the Q-kernel analysis (blue dots) and the virial expansion (black squares), extrapolated to L = oo; The 
Fermi-energy is fixed at the mobility edge. 
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amplitude 

(A(r))2 ^ y 1 d'^vA^r) ^ X J^^ d^T^mliO (167) 

where we used the definition ((TO)) and Eas. (j93l94p to derive the r.h.s. of the above 
equation. Then we calculate its simple average 



A(r) ^ i 1 d'^rA(r) = A dm)^c{0 



(168) 



^From Eas. (jl68ll67p we conclude that 



2 



/=^^^ = ^2(7) = 7^hf «1 (169) 



(A(r))2 C^(.l)\Eo 
where the dunensionless function 



oo 



2 



^^^^"^U f*''''^27;^j ' 2(0-57) -3.1 (170) 

The small ratio (jl69p gives the estimate of the space fraction where pairing correlations 
are well established. Indeed, consider a toy model where A(r) = Ao in a fraction / of the 

whole space, and zero otherwise. Then we find A(r) = /Aq and (A(r))2 = /Ag, leading 
to the ratio (|169p equal to /. In a realistic case A(r) take a continuum of values, so 
one can use this ratio as a proper definition of the fraction of space with well developed 
superconductive gap. 

Note that the regions in space where the superconductive gap is appreciable constitute 
a finite fraction of the entire space, despite the fact that fractal support of any single- 
particle wavefunction occupies vanishing fraction of the entire space. Therefore the global 
spatial pattern of superconductivity is not a fractal but rather is reminiscent of that of 
the multi-fractal metal f7l\ . 

The above estimate was done for a region near Tc- At temperatures much below Tc, 
(|96p should be replaced by the similar equation in which the argument is replaced by 
+ A^, . As one might expect, one ends up with the same estimates in which Tc gets 
replaced by Ao(0). 

Higher moments of the distribution function P[A(r)] can be estimated using the 
algebra of multi- fractal states discussed in section 12.21 see Eq. dTTl) . Straightforward 
generalization of Eq. (jl69p gives 



A(r; 



(A(r))^ 



oc (r,/i;o)<^-''"/'')("-^' (171) 



Therefore the moments of A(r) contain information about all fractal dimensions dn- This 
is because the fractality of single-particle states manifests itself locally at a scale smaller 
than Lt ~ (I'oT'c)"^^'^ , above this length the order parameter and other properties of the 
superconductor become homogeneous, similar to what happens in a multi- fractal metal 



77| which has a sparse fractal structure within the correlation length. Thus, this state 
should be properly named fractal superconductor. 
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Figure 24: Local zero-temperature density of states at the mobility edge of the 3D Anderson model, 
averaged over 5x5 square (blue dots); average T = superconductive DoS at the mobility edge, according 
to Eo. 11731 1. shown by dashed red line; the BCS density of states with the same gap value (full black 
line) . 



5.4- Low-temperature density of states. 

Local density of states (DoS) in the superconducting state is given by 

^(£,r) = IY.(i + ^) [5ie e,) + 5{e + e,)] ^]ir) (172) 



2 ^ V e 

j 

where Sj = + (^j ) and A(^) is the solution of the gap equation (|162p . We begin 
with the calculation of the average DoS: 



de 



(173) 



where the function ^(e) is defined as the (positive- valued) inverse function for e(^) ~ 

We plot in Fig IM] the local DoS (at T ~ 0) obtained numerically using Ea. (|172p . as 
well as the average DoS obtained from Ea. (jl73p for 7 — 0.57; the usual BCS DoS is 
shown for comparison. 

Note that the terms £,j/e in the parenthesis in Ea. (|172p were irrelevant while calculat- 
ing the average DoS because contributions from ^ = ±^(e) cancel out. There is no such 
cancellation for the higher moments of i/(e,r). In particular, these terms might lead to 
an asymmetry of local DoS: z^(e, r) ^ i^(— e, r). To quantify this asymmetry we define the 
antisymmetric part of the DoS: i^a(e,r) = i(z^(e, r) — i^(— e, r). Using then Eas. (|172ll73p . 
we find the variance 
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i^lie^r) = -[uie)^-Y ] [A/(0) - M(2^(£))] (174) 



where the function M(^) is defined in Eq.gl]). Since M(0) ~ > M(2^), the variance 
(jl74p diverges for an infinite system, and one should consider the distribution function 
M^[p(e,r)] for the dimensionless variable 

j/a(£,r)v^e 

The distribution function W^[p(£)] coincides with that of the wavefunction's intensities 
V{ip'^), see Eq.(2.33) of the review [t^- It is determined by the "singularity spectrum" 
/(a); the shape of this function for the 3D Anderson model is similar to the one found 



for the power-law banded matrix model [92|, |72| with the parameter b « 0.4, see Fig. 3 of 
paper [78j. 

It is clear from the definition (|175[) that the normalized asymmetric fluctuations 
i'a(e, r)/:/(e) are small for e « A(0), where ^(e) — )- 0. However, for a generic value 
of the energy of the order of A(0), the asymmetry in the tunneling spectra is strong; 
we emphasize that it occurs due to mesoscopic fluctuations which do not involve any 
"regular" mechanism of particle- hole asymmetry. 



5.5. Superfluid density and critical current. 

Wc define the superfluid density ps via the relation j = —psA./c assuming the trans- 
verse gauge with V • A = 0. To estimate Ps{T) near Tc, one can use the expression (|108p 
for the free energy functional, and make there the replacement V — >■ V — i(2e/hc)A. 
Then taking the double functional derivative S^Fql/SA'^ one obtains: 

Ps{T)^2(^^y i^oT^C^'iT) (176) 

where *(r) = ^l^pj^ = -1^^ = 9.5 (^1 - and the coefficient C - L| is given by 
Eg.dnsl). 

We also estimate the critical current jc in a standard way via the Ginzburg-Landau 
functional ()108p . By extrapolating the result of Ginzburg-Landau theory to lower tem- 
peratures, one obtains 

Jc - c|r,(i.oT,)^/3 (177) 

Extrapolation (|177p of the Ginzburg-Landau result to T = might be wrong by a 
factor of the order of unity, this factor for the fractal superconductor might be different 
from the known factor for a conventional "dirty-limit" superconductor. 

The results Eas. (jl76p . (ll77|) do not contain any "fractal" specificity and coincide with 
estimates of earlier works[^i3| that used scaling arguments. Yet, as we have demonstrated 
above, the web of superconductivity is sparse and occupies only a small fraction / ~ A <C 
1 of space. The derivation of Ea. (|135|) demonstrates that the fractality- independent result 
for the coefficient C emerges because of the compensation of the small parameter by 
the enhancement factor (Eq/Tc)'^'' . We interpret this cancellation as the compensation of 
the small cross-section area of " superconductive filaments" by a large current that they 
support. 
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6. Superconductivity with a pseudogap. 

In this section wo study the superconductivity formed against the background of mod- 
erately locahzcd single-electron eigenstates. We demonstrate here the existence of the 
interesting regime where the local pairing (discussed above in section [3] in relation with 
the hard-gap insulator) coexists with the long-range superconducting correlations. In 
this situation the material demonstrates the " pseudo-gap" phenomenology characterized 
by the insulating behavior of R{T) in a significant temperature range above the super- 
conducting transition. We provide a semi-quantitative description for the evolution of 
the pseudo-gap Ap and Tc while Ep moves across the mobility edge Ec- Surprisingly, we 
find that Tc{Ep) curve is non-monotonic, with a maximum reached at the Fermi-level on 
the localized side of the Anderson transition. The most important outcome of our analy- 
sis is the existence of a range of Fermi energies for which superconductor is characterized 
by relatively high values of Tc and by a larger insulating pseudo-gap Ap ^ Tc- In this 
regime the low energy properties such as the formation of the superconductivity can be 
described by the truncated pseudo-spin model. In this model the effective Hamiltonian is 
of the form (|145p while the Hilbert space contains only the states in which each localized 
orbital is either empty or populated by a pair of electrons. The effective Hamiltonian 
(|145|) does not contain A p because the states that differ in energy by Ap are absent from 
the truncated Hilbert space. These states, however, appear when single particles are ex- 
cited, thus the single particle excitations spectrum contains both correlations and local 
pairing effects. Below we discuss the region of validity of the mean-field treatment of the 
pseudo-spin Hamiltonian (|145p and then turn to analysis of several measurable quanti- 
ties which may provide one with the "proof of the case" for realization of the pseudo-gap 
superconductivity in real material: the temperature evolution of the local single-particle 
density of states i/(e,T), the Andreev conductance spectra Ga{V) at T <C and the 
temperature dependence of the full spectral weight K{T). 

6.1. Transition temperature and insulating gap as functions of Fermi energy. 

We demonstrated in section 15.21 that the transition temperature found by the virial 
expansion to the third order provides results very close to the ones obtained with the 
modified MFA. These calculations were performed for the critical case Ep = Ec where 
the function A/(w) is given by the power law (HH). This coincidence implies that the virial 
expansion to the third order provide the accurate result for the transition temperature. 
Generally, both methods suffer from the neglect of the off-diagonal matrix elements, 
these matrix elements mix the odd and even fermionic sectors which might decrease the 
Tc below the MFA value by a factor of the order of unity (which would be consistent 
with the Gi ~ 1 obtained for the critical case in section B:.2.5p . While there is no general 
reasons why the off-diagonal matrix elements can be completely neglected in the metallic 
and in the critical region, in the region of well localized single-particle states the validity 
of the spin Hamiltonian (|145p is ensured by the large value of the pseudogap Ap. One 
thus expects that pseudospin model provides very good description of the pseudogap state 
and a semi-quantitative description of the crossover from the pseudogap state (where it 
is exact) to the fractal superconductor (where its results might differ by a factor of the 
order of unity). Because the superconductivity suppression deep in the insulating regime 
is due to the decrease in the effective number of neighbors of each pseudospin and this 
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Figure 25: (Color online) Transition temperature as a function of the Fermi-level position, determined 
by the virial expansion method. Blue dots correspond to the pseudo-spin Hamiltonian Eq. ||145| I with 
gy = 0, whereas other data were obtained after including the Hartree-type terms proportional to g|| of 
repulsive (red triangles and black squares) and attractive (green stars and violet diamonds) sign (see 
text for details). All calculations were done at a fixed coupling constant g± = g. 



effect is missing in the modified MFA we shall use the virial expansion applied to the spin 
Hamiltonian ()145p to obtain the properties in both the crossover and pseudogap regime. 

We begin by applying this method to describe the evolution of Tc as Fcrmi-level 
crosses the mobility edge Ec and gets into the region of weakly localized states. In 
this regime, assuming that single electron states are characterized by small participation 
ratio so that corresponding Ap < Tc, we can use the simplified version of the virial 
approach with the distribution function tanh replaced by tanh ^ , as it was done for 
critical states at Ep = Ec in Sec. V. The numerical results for the minimal pseudo-spin 
Hamiltonian Ea. (|145p with g^^ = and coupling matrix Alij generated by random 3D 
Anderson model Eg.ipTj) are shown in Fig 155] by blue points. They correspond to the 
fixed coupling constant g = 1.7 in the original model and, thus, dimensionless coupling 
constant A = gvo varying as function of Ep. The energy dependence of density of states 
is not negligible around Ec as can be seen in FiglH as a result X{Ep) drops considerably 
in the range of Fermi-energies covered in Fig. 1251 Nevertheless, the calculation not only 
confirms our analytical prediction of the enhancement of Tc near the critical disorder 
made in section [4. II but also shows that Tc continues to increase in a significant interval 
of Fermi energies corresponding to the localized region. This surprising result can be 
traced back to the behavior of the correlation function M{uj) = VMij in the domain of 
weakly localized states. Indeed, it was found previously (see the discussion at the end of 
section VIII B of Ref. [t^I) that the effective decrease of fractal dimension ^2, and, thus, 
the increase of 7 = 1 — c?2/c? occurs while going further into the insulating domain. For 
small values of A it results in the increase of Tc oc A^/t evident in FigHHin the range of 
Ec < E < 6.5. Note that in contrast to Fig[23]the Fig|25] shows the data for the fixed 
value of system size, L = 20, with no finite-size scaling adjustment. This is why the 
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values of Tc data for the critical states, E = 5.5, differ slightly between these two data 
sets. 

The generic pseudo-spin Hamiltonian Eq. (|145| ) contains also the Hartree terms 

-ffHartreo — ^ .9|| MjjS^S^. (178) 

ij 

Moreover, for the short-range 6-function two-body interaction, the coupling constant 
is equal to the Cooper-channel coupling constant g± and, thus, the terms containing the 
spin-spin interaction in pseudo-spin Hamiltonian (jl45p become isotropic. Because the 
full Hamiltonian (jl45p remains highly anisotropic due to the presence of the single spin 
terms, the effect of the S~ and S^Sj spin-spin interactions is dramatically different: 
the former leads to the superconducting instability at low temperatures, while the latter 
term by itself has only weak effect of spin susceptibility for small g^^iy. Because the 
transition temperature depends strongly on the spin susceptibility, in the presence of the 
S^S~ the Hatree term can have slightly more significant effect. Below we study the 
change in the critical temperatures induced by the SfSj interaction for a generic ratio 
= g\\/g± using the virial expansion method. 

The results are shown in Fig.[25]by green and violet points for a = -1-0.5 and a = +1, 
and by red and black points for a = —1 and a ~ —0.5, respectively. We conclude that 
upon the account for the Hartree terms, all the qualitative features of the original Tc{Ep) 
behavior, including the enhancement of Tc by disorder, are preserved. Therefore below 
we stick to the simplest version of our model with a — Q. 

With a further shift of Ep into the localized domain, the two effects become impor- 
tant. First of them is the growth of the inverse participation ratio Mj, and the related 
development of the parity gap Ap , see Ea. ([T7|) . In some range of Ep, as we will see 
soon, the typical value of Ap becomes larger than Tc, which means that the calcula- 
tions of the virial coefficients should be modified, taking into account exactly the effect of 
" dilution" by the odd states discussed in section 14.41 The necessary modification of the 
formalism is presented in the [Appendix B[ The dilution elimination by the pseudo-gap 
leads to the increase of Tc values (as compared to the data for Ap = shown in FigESl 
where the maximum effect of dilution was assumed). However, as one can see in Fig. 
in the region where superconductivity exists this increase is not more than 30%. 

The second effect, which takes over sufficiently far in the localized region, is due to 
the drop of the effective coordination number Z^s of the spin model defined by Eq. (|145p . 
Indeed, the effective number of states k that are coupled to a given state j can be roughly 
estimated by Z ~ voTcLf^^ (this estimate misses the important logarithmic factor that 
we discuss below around Ea. (jl79|) ). With decreasing of the localization length Ljoc, 
the coordination number Zcs drops down, and eventually becomes less than unity, which 
makes description in terms of average matrix elements VMij = M{u}) meaningless. When 
the Fermi energy is increased further, Tc starts to drop down sharply simply because most 
individual pseudospins (with a possible exception of rare small spin clusters) become 
essentially decoupled from each other. 

Evolution of both Tc and A p with the Fermi- level position are shown in Fig|26l where 
we show the results obtained by the virial expansion method, with a finite pseudo-gap 
Ap taken into account. The most important feature seen in this plot is the existence 
of some range of Ep where superconductivity with an appreciable Tc exists, but the 
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Figure 26: (Color online) Virial expansion results for Tc{Ep) (red squares) and typical pscudogap A 
(black dots) as functions of Ep . The model with fixed value of the attraction coupling constant g = I 
was used; pairing susceptibilities were calculated using equations derived in Appendix B. 




Figure 27: (Color online) Virial expansion results for Tc(Ep) (red dots), the typical pseudo-gap A 
(black squares) and the corresponding level spacing Sj^ (green diamonds), as functions of Ep in tl 
semi-logarithmic scale. 
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parity gap Ap exceeds substantially the superconductive gap A « 1.7Tc. Below we 
concentrate specifically on this range of parameters (approximately, it corresponds to 
the range Ep £ (7.5,8.5) for the 3D Anderson model studied here), where the specific 
features of superconductivity developed on top of the pseudo-gap are most pronounced. 
This is precisely the region of the pseudo-gap superconductivity where the pseudo-spin 
model Ea. (|145D provides a quantitatively correct description of superconductivity. 

Note that in this range of parameters the typical level spacing Sl inside the local- 
ization volume (estimated with the use of Eg. ([71]) ') is much larger than both Tc and Ap, 
as it is demonstrated in Fig |27l The very existence of a nonzero Tc in this situation is 
unexpected. [f| It is related to the fact (noticed in Ref. 77 [ and discussed in section !^. 2. 5p 
that the localized states characterized by the localization length Lioc overlap effectively 
up to the distance i?(w) ^ Lioc between them, where a; <C (Sl is the relevant energy 
window (—Tc,Tc). In the present context we can put u ^ Tc and estimate the effec- 
tive interaction range Rq and coordination number Zcs of the spin model (I145D using 
the definition in Eq. (p5)) in subsection 12.2.51 see also Eqs.(21,25) in the paper [77[. We 
obtain 

(179) 



47r g 

ZcS = —Rq{Tc) ■ 2Tc 



where the function Ro{uj) should be determined numerically. Using data shown in Fig. 1271 
for the case of Ep = 8 and similar data for several other values of Ep in the pseudogap 
range, where Tc <C (5l , we obtain the following values for the effective number of neighbors 
in the energy window (—TcTc): 



Ep 


7.75 


8.0 


8.25 


8.5 


Tc 


0.057 


0.044 


0.030 


0.02 


Sl 


1.0 


1.9 


3.5 


6.5 


ZcS 


0.08 


0.09 


0.09 


0.07 



(180) 



Of course, the numerical coefficient in Eq. (|179p is somewhat arbitrary; the main con- 
clusion from Ea. (|179p is that (due to a presence of the logarithmic enhancement in 
i?o(w)) the effective number of neighbors remains nearly constant in a wide range of 
small Tc/5l ^ 1. Together with a weak Tc/Sl dependence this explains why a signifi- 
cant Tc is found even for 6l/Tc ~ 300, see Fig. 

As we explain in more detail below, the possibility to find a pseudogaped super- 
conductive state with a numerically small effective coordination number {Zcs ~ 0.1) is 
natural. The qualitative reason is that our definition of ZcS implies average over all 
states. Therefore, this average includes the states that overlap strongly with each other 
and the ones which overlap little with any other states. The former states form an infinite 
superconducting cluster while the latter do not contribute to superconductivity and are 
largerly irrelevant. In the strongly disordered regime the superconducting cluster occu- 
pies a small part of all states, so the average is dominated by the states which overlap 
very little with each other. 

Upon a further increase of the Fermi energy into the localized region the effective 
number of neighbors slowly decreases leading to a sharp drop in T^yP{Ep), see Fig. 



"^In contrast, superconducting correlations in ultrasmall metall grain are destroyed at 5 ~ A, see 
discussion in section [2. 1.1 1 
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Finally, when the effective number of neighbors becomes less than some critical value 
Zc, the superconducting instability ceases to occur even at T = 0. The details of this 
quantum transition and the qualitative properties of the phase formed at high disorder 
are discussed in [g^ . This work solves the simplified model similar to (jl45p in which spin 
arc located of a Bcthe lattice with coordination number ^Bcthc^ -^^ ^^^^ model Afy — 
j^^^Bctho jTQ^. ^Botho neighbors and otherwise. The result of the solution [93] which 
is important for the present discussion is that critical fluctuations become important 
only in the narrow vicinity of the quantum critical point at which superconductivity 
disappears completely. In other words, the critical fluctuations become important only 
in the systems in which transition temperature is already strongly suppressed by disorder. 



In more detail, the work |93| distinguishes between the high temperature phase where 
corrections to mean field are small and a low temperature phase where mean field solution 
becomes qualitatively incorrect due to the strong inhomogeneity of the formed ordered 
phase. The temperature (Trsb) that separates them corresponds to replica symmetry 
breaking (RSB) in the formalism, it turns out to be numerically small even for modest 
number of neighbors as illustrated in Fig. Viewed differently, it means that the 

effective number of neighbors in this regime is very low (see Fig. 1281 insert). 

As long as the transition temperature obtained in the mean field approximation is 
above Trsb shown in Fig. 1281 the mean field equations are qualitatively correct. Thus, 
one can use mean-field approximation for the Hamiltonian (jl45p to get semi-quantitative 
results for the transition temperature in a large interval of Fermi energies at which Tc{Ep) 
decreases significantly from its maximum. 

The solution of the Bethe lattice model [9^ can be used to estimate critical value 
of the effective coordination number = vbTcq ^rsb'" which corresponds to the 

value ^^gB^° of the Bethc-lattice coordination number where Tc = Tbsb and standard 
mean-field solution becomes qualitatively incorrect. The details of the corresponding 
calculation will be presented in the extension [9^ of the paper [93|; the result, valid in 
the limit Zecthc S> 1, is 



^Bethc X ^l/2\b,, rp 1 rZ„-C/2 Vi^bTrSB 

^eff =A.e/ .BTeo = -V^e |i,g(,^j.^^^)|«l (181) 

where A;, = z/b^b ^ 1 is the dimensionless coupling constant for the Bethe lattice model, 
with vbTcq ~ ■|e*^e~^/^'', and C = 0.577.. is the Eulcr number. Ea. (jl8ip demonstrates 
again, that typical number of neighbors with energies within the energy stripe ^ can 
be very small in the weak-coupling limit, contrary to naive expectations. Physically, the 
reason for this is that in this regime a dilute superconducting cluster is formed by a 
small number of sites. One expects that similar phenomena should happen in physical 
systems. Although our fractal model differs significantly from the Bethe lattice model due 
to power-law dependence of the interaction strength M[uj), we expect that qualitative 
conclusion about 2'^'" < 1 survives in the fractal model as well. 

The significant number of neighbors that persists up to numerically large values of 
5l/Tc ^ 300 implies that Ginzburg parameter remains of the order of unity in this 



'''Notice that this definition of the number of neighbors differs a lot from the II179II that counts the 
neighbors only in the narrow energy interval of Tc- 
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Figure 28: Temperature separating high and low temperature regimes of the simplified pseudospin model 
11145 II on Bethe lattiee with ^sethe neighbors, for which Mij = l/^Bethe a-iid density of states per site 
i/g. In the high temperature regime T > Tj^sB the mean field theory is semi-quantitatively eorrect, 
in the low temperature regime it is qualitatively wrong; the value of Tugg does not depend upon 
coupling constant g within the Bethe lattice model. The inset shows the effective number of neighbors 
in this model defined by Z^"^"^ = J^BTco^Bethe which the mean field solution becomes qualitatively 
incorrect. 



range. In this regime the gap below Tc can be semi-quantitatively approximated by 
the solution of the mean field equation: 



where we used the fact Ap ^ Tc to replace 2T T and M(aj) is given, for localized 
states, by interpolating formula, Ea. (|5n|) . 

One of basic results of weak-coupling BCS theory is the universal relation between 
low-tempcrature gap value and transition temperature, A(0) — l.TQTc- As was shown in 
section [5TT1 the ratio A(0)/Tc stays about the same within MFA for fractal superconduc- 
tor with 7 exponent near 0.6, corresponding to 3D Anderson model, see Figl^Db. The 
presence of strong pseudo-gap Ap ^ Tc leads to the doubling of Tc for the same value of 
A(0); on the other hand, thermal fluctuations beyond MFA always lead to suppression of 
Tc with respect to A(0). We do not attempt here an explicit calculation of the A(0)/Tc 
ratio in the pseudo-gap region, which should take into account both the above effects; 
roughly we may expect A(0)/rc ^ 1-5 — 2 in this region. 



*This does not cotradict our findings in section [4.2.51 that Ginzburg parameter becomes Gi ~ 0(1) 
already for the fractal superconductor with the Fermi energy at the mobility edge. It only implies that Gi 
grows very slowly and remains 0(1) upon a further increase of the Fermi energy untill the the transition 
temperature becomes small. The reason for it is the appearance of a large pseudogap Ap S> Tc that 
suppresses the fiuctuational processes that involve off-diagonal terms Mijf^i. 
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Figure 29: Distribution functions P{Tc) obtained by the virial expansion for different Fermi energies. 
All values of Tc were obtained in the pseudo-spin approximation Ap ^ Tc, which leads to some overes- 
timation of Tc's at the Fermi energy Ep = 7.9 where Ap/Tc is not very large. The bi-modal character 
of the distribution for 8.2 < Ep < 8.7 indicates on the percolation character of the superconductive 
transition. 

With a further shift of the Fermi energy deeper into the locahzed region {Ep > 8.2 for 
our model), a remarkable new feature appears in the distribution of "transition temper- 
atures" obtained from the truncated criterion Eq. (|150p via the virial expansion. Namely, 
it acquires two maxima as shown in Fig. 1291 One of them is located at Tc = and corre- 
sponds to an insulating phase, while another one corresponds to a superconducting phase 
with a nonzero Tc- The magnitude of this latter maximum decreases with increasing Ep, 
until it disappears completely at about Ep = 8.7. 

Dramatic modification of the character of the distribution of Tc implies that the 
qualitative change in the physics of the superconductive transition. The appearance of 
two maxima in the distribution function means that not-too-large systems available for 
numerical studies are either superconducting or insulating. Larger systems and physical 
materials become strongly inhomogeneous with some parts becoming superconducting 
and some insulating at low temperatures. This observation is fully consistent with the 
result of the analytical theory[93| that predicts a replica symmetry breaking transition 
upon the decrease of the effective number of neighbors. This transition signals the absence 
of a local self- averaging and a global behavior which is dominated by rare superconducting 
paths. Thus, it is likely that for the physical, very large system the transition in this 
energy interval occurs by a percolation scenario. 

The inhomogeneity of the superconductor in this regime is very different from that of 
the "fractal" superconductivity that we discussed in the previous section. The latter is 
macro scopically homogeneous because the strong inhomogeneity corresponding the frac- 
tal structure is seen only at scales smaller than Lt- In this case, the distribution function 
of transition temperatures for all macroscopic samples with L > Lt appears to be the 
Gaussian single-peak distribution with the width decreasing to zero as the system size 
increases, typical of the self-averaging quantities. In contrast the bimodal distribution 
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observed in numerical simulations and in analytic solutiont93| of the pseudospin model 
signals of the macroscopic inhomogeneity and lack of self averaging. 

Finally wc note that our treatment above assumed a second-order transition to the 
superconducting phase, for instance, wc have determined the transition temperature as 
a temperature of Cooper instability. We cannot exclude the first order transition to the 
superconducting phase although it seems unlikely, especially in the regime of a large 
disorder. 



6.2. Tunneling conductance. 

Here we analyze in detail the tunneling conductance into the pseudo-gaped state, at 
temperatures above and below the superconductive transition. We compute the average 
tunneling conductance in these two cases in sections I6.2.1l and l6.2.2l These conductances 
can be probed by large tunneling junctions. Computation of local tunneling conductance 
(as measured by scanning tunneling probe) is more complicated, we are able to make 
only qualitative predictions for this quantity in section [6.2.31 

6.2.1. Tunneling in a normal state. 

We start by discussing the temperature-dependent differential tunneling conductance 
G(y, T, r) in the presence of a pseudo-gap Ap above the superconducting transition. In 
this case the total tunneling current is given by 

I{V) = / /dedeiz.o(e,r)[(l-/(£))/(ei)5(si-e + ey-A«) 

evo J J L 

-(1 - f{e,))f{e)5{e, - e + eV + A*^')] (183) 

where Go is the bare conductance of the tunnel junction, i/o(£,r) = "Y^- tjjj{r)6{e — ^j), 

andlj =0-A^'\ Eg.lfTMl) was obtained by the following arguments. First, we consider 
the tunneling from the probe to the sample, it is given by the first term in (151). Here 
the j-th level in the sample is empty, the probe level is full, which gives rise to the usual 
factor (1 — f{e))f{ei) involving the Fermi distribution function /(e). In this case the 

energy conservation requires £i = £,j — eV = e + Ap^ — eV. Second, we consider the 
tunneling from the sample to the probe: the j-th level in the sample is full, the probe's 
level is empty, the thermal distribution functions produce the factor /(e)(l — /(ei)); in 
this case the energy of a pair of electrons on the j-th level contains the binding energy 
~2Ap\ so that the energy conservation reads: £i = — 2Ap'' — eV ~ s — Ap^ — eV. 

Note that this calculation assumes the validity of Fermi Golden Rule for the tunneling 
rate; it also neglects electron-electron interactions (apart from the local terms leading 
to the parity gap Ap). The simultaneous validity of both these assumptions is not 
guaranteed when we consider localized states with relatively large typical level spacing 
5l', we will discuss this more below in section [6.2.31 

The differential conductance G{V) — dl/dV corresponding to Eg. pSSp can be repre- 
sented in the usual form 

G{V,T,r) _i , ,f df{e-eVy 



Go 
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where the density of states in the presence of a pseudo-gap (but above Tc) depends 
expHcitly on temperature: 

Me,r) = ^vl(r) [5(.-|,-A(^')(1-/(£-A«')) (185) 
j 

+ 5(e-|,+A(,^')/(e + A^')] 

Note that Ea. (|185p cannot be represented in terms of the normal-state DoS, i^o{e, r) , 
due to the correlations between the eigenfunction intensities 7/'|(r) and the local pairing 

gaps Ap''. The spatial average of the tunneling DoS can be obtained from Ea. (|185p by 
integration over r which replaces of Vip'j{r) by unity in Ea. (|185p . 

The next (approximate) step in the simplification of this expression is to replace the 
summation over eigenstates j by averaging over distribution of pairing gaps P{Ap)dAp, 
see Fig. [HI The distribution function of the normalized gaps y = Ap/Ap"^ can be fit by 
analytical expression 

P{y) = ^cxp(---hy-b2yA (186) 

where c = 1.1, hi ~ 0.16, 62 = 0.03 and A = 1.8 for E = 8.0 (coefficients 61,62 are not 
universal as can be seen from Fig. [T5)) . 

The result for average density of states then reads: 

cosh,. + cosh f ^(^^-^^ 

where z = Ap^/T and i'o(e) is the average (in general, energy-dependent) DoS in the 
normal state. In the low-temperature limit Ea. (jl87p reduces to Eg. ([75)1 from section [3] 

Evolution of the average tunneling DoS obtained from Ea. (|185p and its approximation 
Ea. (|187p at different temperatures T < A^p^ for the 3D Gaussian Anderson model are 
shown in Fig l30l The asymmetry with respect to the sign of e is due to a non-negligible 
energy- dependence of the bare DoS i'o(£). 

The obtained results for the DoS are translated into the measurable tunneling con- 
ductance by Ea. (|184[) . The zero-bias conductances G(0,T) obtained from the "exact" 
fEa. ()185[) ) and the approximate (Ea. (|187p ) expressions for the DoS are shown in Fig.[3T] 
It is seen from Figsl30l and [311 that the approximation ()187p works reasonably well; below 
we will use this approximation in the further analysis of the average tunneling conduc- 
tance. 

The ensemble-averaged tunneling conductance G{V, T) as a function of voltage is 
plotted in Fig. [32] for several temperatures. Remarkably, the curves we obtained demon- 
strate nearly exact crossing at eV ~ Ap^. In other terms, the tunneling conductance at 
this voltage is nearly T-independent. This unexpected feature provides a simple way for 
experimental determination of Ap"^. Note that at very low temperatures T <^ A^p^ the 
crossing point moves to somewhat lower voltages. Indeed, a simple analysis of the integral 
in Ea. (jl87p shows that in the T — > limit the derivative dv{e)/dT oc dP{£)/d£ vanishes 
at £ = Eo corresponding to the maximum of the function P(£/Ap''), i.e. £0 ~ O.SAp^, 
see Ea. (fT86l) and Fig. [H 

The averaged tunneling conductance can be measured in the normal state by a large- 
area tunnel junction contact. It is important that the temperatures are not too low 
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Figure 30; (Color online) Global tunneling density of states for the 3D Gaussian Anderson model at 
Ep = 8 and the system size 20'^ for different temperatures normalized to the typical value of the local 
gap Ap^. The data points stand for the all- numerical evaluation of Eq. lll85ll whereas the lines were 
calculated by means of Eas. lll87ll and II186I I. 
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Figure 31: (Color online) Zero-bias tunneling conductance as a function of temperature below the 
typical pseudo-gap value; results obtained with exact formula {TSS} are shown by black dots, whereas 
blue squares correspond to the approximate formula II187II . 
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Figure 32: (Color online) Ensemble- averaged tunneling conductance, Ea. lll84ll . for several values of 
temperature below the typical local gap Ap^. 



with respect to the typical Ap because as temperature drops down, the bulk resistivity 
grows exponentially, leading to development of a strong Coulomb gap (not taken into 
the account above) which makes tunneling measurement unfeasible. The same problem 
becomes even more severe for the STM measurements at low temperatures. 



6.2.2. Tunneling in a superconductor. 

Provided that Fermi Golden Rule remains valid the tunneling current in the super- 
conductive state is given by the formula similar to Ea. ()183p . We only need to replace 
in Ea. ()183p the bare normal-state DoS i'o(e, r) by its superconductive counterpart given 
by Eq. (|172p . Then, repeating the derivation of Eq. (|185p . we come to the expression for 
the tunneling density of states in presence of both pseudo-gap and the superconducting 
correlations: 



1 




x(l-/(£-A^"))) + 



A 



<5(e 



die - Sj + A^'^) + 6ie + e, + A'^>) f{e + A'^>) 



0> 



(188) 



where 



(II82|) 



J J and the function A(^) has to be determined by solving the 
Below we use Ea. (|188|) 



gap equation (jl82|) . Below we use Ea. (|188|) in order to derive an expression for the 
averaged tunneling conductance; situation with the local tunneling conductance is more 
complicated and cannot be described this way, as explained in subsection 16.2.31 
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Figure 33: (Color online) Ensemble- averaged tunneling DoS, Ea.l ll89l l. for several values of temperature 
mueh below the value of the superconductive gap Aq = Ao{5 = 0). The latter is chosen to be Aq = 
O.SAp'', where Ap^ is the parity gap. 

The averaged density of states in the superconductive state with a pseudo-gap, h'sds), 
can be obtained in the same way as it was done above for the normal state, see Eq. (jl87p . 
The odd terms in proportional to ^j/(e± Ap) in Ea. (|188p cancel out upon summation 
over S,j = provided wc neglected the energy-dependence of typical Ap(E) inside 

the energy band of the width AE ^ Ap. We get 

/•OO 

iysc{e)= / P{y) 
Jo 

where z = Ap^/T and i^sc''(e) is the DoS in the superconductive state without a pseudo- 
gap, defined by i/i?(£) = iyo{s)\d^{e)/de\ (|173p . similar to the one shown in Fig. [Ml for 
the case of critical superconductor. Note that the function ^(e) that we need here differs 
slightly from the one shown in Fig. l24l due to the different form of the correlation function 
M{u!) entering modified MFA equation (|182p in the case of pairing between localized 
states. However, the effect of this difference reduces to an overall prcfactor, as is seen 
from Fig. [51 We thus do not expect significant difference between the solution A(^) and 
the corresponding ^(e) functions calculated within the modified MFA approximation in 
the mildly localized region considered here and for the "critical" case studied in section 
m Thus we use d£,{e)/de\ as it was found in section [5^ but take into account a smooth 
energy-dependent DoS pre- factor ^o(£)i as it was done above when computing Ea. (|187p 
and Fig. [201 

Evolution of the z^sc(e) shapes as function of temperature in the low-temperature 
range T < O.STc is shown in Fig l33l where we use the analytic interpolation formula 
Eg. pSOp for the distribution function P{y), and the T = solution for the gap function 
A(^). The average tunneling conductance G{V, T) obtained with this DoS from Ea. (|184p 
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Figure 34: (Color online) Average tunneling eonductance in the superconduetive state for several values 
of temperature much belowr the value of the superconduetive gap Aq = 0.5Ap^. 



is shown in Fig. 1341 The crossing point at this plot, i.e. the T-independent conductance 
at some specific vahic of voltage, corresponds to eV — Aq ~ 0.5 A^p^. In zero energy limit 
the expression for the density of states simplifies to 



ev/T + 1 



dy 



(190) 



The temperature dependence of the corresponding zero bias tunneling conductance G(0, T) 
is shown in Fig ]35l 

The averaged DoS and the tunneling conductance presented in Figs. [SH |35] do not 
show any coherence peaks similar to the one presented in Fig. [53] for the critical su- 
perconductive state (section 15. 4p . This is due to the averaging of the peak position 
eV = Aq + Ap over a broad distribution P{Ap) of the local pairing gaps A^'^ Note 
that the non-monotonic behavior of the DoS curve in Fig. [32] reflects the (smeared) 
edges at e = Ap^ ± Aq of the i^s? (e) function shifted by the paring gap Ap ^ > Aq and 
weighted by the temperature-dependent distribution functions in Ea. (|189p . 

We emphasize again that the results obtained in this section give the averaged tun- 
neling conductance which can be measured with a large-area junctions, such as the ones 
used in |4l| . They may differ significantly from the local tunneling conductance measured 



by a small tip and discussed in the following subsection. 
6.2.3. Point contact tunneling. 



The equation (|183p can be easily generalized to describe point conductance. One 
gets using the perturbation theory in the tunneling amplitude, t (which is equivalent 
to Fermi Golden Rule) that the most general case the point conductance is given by 
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Figure 35: Zero-biased conductance as function of temperature for the case of Aq : 
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equation p83p where i'o{e,r) is replaced by the exact single particle Green function at 
the tunneling position, IniGfl,(e,r). In a poor conductor in which single particle states 
are localized or almost localized one expects that ImGi?(e,r) shows a series of peaks 
corresponding to the single particle states located in the vicinity of the tunneling point. 
Applying this reasoning to the superconductor-insulator transition in a situation when 
Ep corresponds to a well developed paring gap Ap and repeating the arguments of the 
previous subsection one expects to observe the peaks in the tunneling conductance with 
width T and separation Sl- 



These arguments do not take into account the effect of the collective modes of the 
superconductor which smear the peaks in the tunneling conductance. In fact, the pres- 
ence of the collective modes make the observation of the peak structures in the density 
of states difficult. Briefly, the experimental observation of tunneling is possible only in 
the superconducting state, but in this state the effect of the collective modes on the sin- 
gle particle density of states is large. We estimate this effect within a simplified model, 
namely we use the pscudospin Hamiltonian (jl45p to which we add the possibility of 
having single particle excitation at the single site 0; 



H = ^on + Ap{l 



j 



HpsiS) (191) 



where n = J2a '^a'^a — 1 where a =t, i are two spin components, and the sum over j runs 
over the neighbors of the site 0. This model describes the process of electron tunneling 
in a single localized state; in a more realistic case the electrons tunnel in few such states 
with an amplitude proportional to ipai''~) where r is a tunneling point. These processes 
are incoherent due to randomness of local binding energies Ap'' , so the physical tunneling 
conductivity is the sum of the contributions from different eigenstates. 
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In a good BCS-like superconductor characterized by a large number of neighbors one 
can replace pairing field h± = g J^j ^ojSf by their average values {h±_) = A. In this 
case the Green function of the electron at site obeys Dyson equation 

GrU) = ^- (192) 

Gm{e) = ,^-n - (l^^) 

e-Co- Apsgn(Co) + iO 



The resulting single electron Green function has a pole at eo = \/ (Co + ^psgn(^o))^ + 
which shows the combined effect of the pseudogap and superconducting order parameter 
on the electron. As one expects, in the case of a large number of neighbors the electron 
density of states has a sharp peak at the pole energy eo which is smeared only by thermal 
effects by 5e ~ T . As the disorder is increased the effective number of neighbors goes 
down and fluctuations of the pairing field become relevant. This leads to two physical 
effects. The first effect is that tunneling of the electrons removes the spin from site 
which changes the effective spin Hamiltonian. The new ground state spin wave function 
is orthogonal to the one realized before the tunneling process. This leads to a suppression 
of the tunneling of the electron with the energy exactly equal to eq- The second effect is 
that tunneling might be accompanied by the emission of the pseudospin excitations. This 
leads to non-zero tunneling density of states at energies larger than eo. Both effects are 
controlled by the same parameter, the effective number of neighbors which determines 
the strength of the quantum fluctuations of the pairing field and both round the peaks 
in the observed spectrum. In order to compute local tunneling conductance one needs 
the full description of the spin fluctuations in the regime of the small Z^f / which is not 
presently available. 

Below we estimate the second effect, i.e. the magnitude of the level broadening in 
the framework of the perturbation theory. In perturbation theory the emission of the 
pseudospin excitations is described by the inclusion of the self energy part in the Dyson 
equation (|192p : 

Gj{{e) - A2GA(-e) e - Co - Apsgn(Co) - Sj?(e) 



ImEfl(e) = .g'^Moj-y^ -^linG a{<^ ^ e)lxnD r{u 



(195) 



where D{uj) = ((5^5^)) is irreducible correlator of the pairing field. Even in a weakly 
superconducting phase the pseudospin excitations are delocalized and thus have a con- 
tinuos spectrum. Here we focus on the large energy scale E ^ Ap properties of the 
spectrum, we do not discuss here the structure of tunneling conductance very close to 
the edge, where coherence peaks are expected. At these large energies, the density of 
collective modes is featureless, so we estimate the lmD{uj) = VBd{uj). 

As it is clear from (jl94p . the interaction with pseudospins lead to the decay rate of 
the single particle excitations, Vf,i , which smears the peak structure of ImG(e) replacing 
it with the threshold behavior. In the leading order of perturbation theory we replace 
G(a; - e) in pM| by Go(w - e) and get for the level width T^i = g'^VB J2j ^oj- 
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Because the expression for the level width, T^i contains the square of the interaction 
constant, g^, but only one sum over neighboring sites, it is small in a good supercon- 
ductor (where gM cx 1/Z) but may become large in a weak superconductor close to 
superconductor-insulator transition. We estimate it in the framework of the Bcthc lat- 
tice model for this transition discussed in section [6.11 In this model Moj = 1/ZBethc, so 
we get Fei = g'^i^B /ZBcthc. lu order to express it in terms of the physical parameters, we 
use the equation (|18ip that relates the transition temperature to the effective number 
of neighbors 2^Bcthc. at the onset of the weak superconductor regime. The onset of this 
regime can be defined (and determined experimentally) by the beginning of the sharp 
decrease of Tc with the increase of disorder. We get: 



Tel - {gi^BfVT^^ T, . (196) 

The apparently large value of the level width in this regime implies that distinct peaks 
become absent, being replaced by the threshold. For a tunneling leading to a single 
state, a = 0, the threshold coincides with the value of energy « [^^ -I- A^] ^^"^ -\- A^p\ 
This discussion assumed that tunneling happens only in one localized state; in a more 
realistic situation when tunneling amplitudes proportional to ■(/'q(^) are significant for a 
few states, the threshold energy is equal to A plus a minimal value of min Ap available 
for a contact tunneling process at a given point. 

We conclude that the model Hamiltonian Eq. (|19ip results in the threshold depen- 
dence of the tunneling conductance. The threshold distribution (observed when scan- 
ning along the surface) is expected to be similar to P(Ap) distribution shown in Fig. 1141 
whereas the tunneling conductance G{y) should be similar to the integral of that dis- 
tribution shown in Fi g. | 1 5l These conclusions are in a qualitative agreements with the 
data presented in (i^ lQSf. 

The model p91D neglects the off-diagonal matrix elements of the form 



that are present in the full Hamiltonian and might be important for the quantitative 
description of the data. Physically these terms lead to the possibility of non-local process 
in which electron is converted into a hole in the adjacent state and emits a collective mode 
quantum. The presence of a sum over many states j connected with a pre-selected state 
a leads to a partial average over the gap distribution P(Ap); this suppresses the relative 
magnitude of spatial fluctuations of the local thresholds. As a result, we expect that the 
threshold voltage is given by 

eT4h = athAp + A (197) 

with ath < 1 being numerical factor of order unity, with a moderate spatial fluctuations. 

The arguments leading to an estimate Eq. (jl96p imply that peak structure of the 
density of states might be only observed in a narrow transient regime that corresponds 
to a well developed pseudogap combined with a strong superconductor. It is not clear 
that this re gim e can be realized in physical systems. In particular, it is very likely that 
experiment 42, 9^ which reported the gap in the density of states but no obvious peak 



structure were performed in a weakly superconducting state close to the transition into 
the insulator. As explained above, in this regime one expects that at a typical point 
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Figure 36: (Color online) Andreev point-contact conductance (normalized by its value at zero bias) 
is calculated within zero-temperature BTK approximation, Ea. lll98ll . for three different values of the 
contact transparency, 7" = 0.2 (blue line) , 0.4 (red) and 0.6 (yellow). Voltage eV is normalized by 
the maximum value Aq of the gap function A(5). It is assumed that single-electron tunnelling is fully 
suppressed by large gap Ap 2> Aq. 

the spectrum of collective modes is continuos and featureless resulting in the tunneling 
conductance characterized by a threshold that corresponds to the minimal value of eq. 
However, large inhomogeneity of the ordered state (discussed in a different context in 
the section 16. ip implies that at some locations the spectrum of collective modes might 
develop large peaks leading to the peak structures in the tunneling conductance. 

Qualitatively, interaction of single electron with pseudospin collective modes plays 
the same role as its interaction with other electrons in zero bias anomaly and Fermi- 
edge singularity problem The quantitative theory of these effects and in particular, 
the predictions for the shape of the pseudogap require the full theory of the collective 
modes appearing in the vicinity of superconductor-insulator transition, this is beyond 
the approach developed in this paper. 

6.3. Andreev contact conductance at low temperatures. 

In the regime of a strong pseudo-gap Ap > Ao a peculiar situation occurs: the 
energy gap for a single-particle excitation Ai = Ap + Aq is larger than the two-particle 
excitation gap 2Ao. The crucial point here is that the local paring energy Ap drops 
out of the expression for the two-particle excitation gap which is determined only by 
the collective superconductive gap Ag. Local paring energy becomes important only 
when we have to break the pair to produce the single-particle excitation. This allows 
to experimentally distinguish between the two types of gaps: the paring gap Ai seen 
in the single- electron tunneling and the collective gap Aq which can be observed in the 
experiments where the entire pair is transferred through the tunnel junction. This type of 
tunneling experiments is called Andreev-contact spectroscopy (see (97| for the review of 
its application to similar studies in high-Tc cuprates). Since the probability to transfer a 
pair is proportional to the square of the tunnel contact transparency T, the latter should 
be not too small for the Andreev conductance to be observable. 

Below we provide simplest estimates of Andreev conductance for the pseudogapped 
fractal superconductor, based upon Blonder- Tinkham-Klapwijk (BTK) approach [98|. 
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We consider the case of moderately large parameter of normal reflection Z > 1. Charge 
transport within BTK approach is characterized by (energy-dependent) probabilities A, 
B, C and D of four processes which may occur when normal incident electron approaches 
N-S boundary: A{e) is the probability of Andrccv reflection, i?(e) corresponds to nor- 
mal reflection, C(e) and D{e) correspond to two different processes of single-electron 
transmission into superconductor. We are interested here in the low-temperature limit 
T "C Tc ^ Ap and relatively low voltages eV < Ap, thus the processes of electron 
transmission into superconductor can be neglected and the full sub-pseudogap current 
can be written as 

Ins{V) ^ I de 2A(e)[/(e - eV) - /(e)] 



This expression is obtained from Eq.(17) of the BTK paper by the substitution of the 
integrand \ + A — B = 2 A + C + D ioi 2 A, due to vanishing of normal current. There- 
fore differential conductance is proportional to the Andreev reflection probability A(e), 
which we will write using BTK results, see Table II of Ref. [o^, but in slightly different 
notations. To begin with, we introduce normal transmission coefficient T = 1/(1 + Z'^) 
and write NS differential conductance in BTK approximation (at T = 0) for usual su- 
perconductor as 

^ = G^siV) (198) 
■AieV)Y nlieV) 



eV J T'nlieV) + {2 - T)^ +2noieV)i2- T)9{eV - A) 



where Gg = 2TNch is the normal-state conductance of the contact, Nch is the number 
of orbital conductance channels (in BTK approach all channels are characterized by the 
same transmission coefficient T), and no{e) = e/\e^ — A^|^/^. 

In order to generalize Ea. ()198p for the case of fractal superconductor, we consider 
subgap region e < Aq and region of real excitations e > Aq separately. In the subgap 
region we introduce function (e) defined by 

--(^^-U'^ -e^ + e + AHO ^'''^ 

where A(^) is the solution to the gap equation (|182p . Then NS conductance in the subgap 
region is: 

G-^siv) = 2Gor^^2_(,^)„i)-^2 + (2_r)2 

For the region above the gap one should introduce function A(e) = A(^(£)), and replace 
in Eq. (jl98p the function no(e) by the actual normalized density of states n^{e) = v{£)/vq^ 
defined by Ea. ([T75)) : 



r+ m 2C r {HeV)/eVf nl{eV) 



The comment is in order here: wc do not expect Eas. (|200l20ip to be quantitatively exact 
for the problem considered, due to a number of oversimplifying assumptions borrowed 

90 



from the BTK approach (it was developed for clean superconductors in a contact with 
clean metal). However, comparison of BTK expression for A{e) with exact calculations 
done for dirty NS ju nctio ns (see Appendix to the paper [o^ and references in that pa- 
per, and also paper [lOOf ) shows that the key feature of the result (|200l20ip . that is, 
the peak of differential conductance at eV = Aq, is reproduced by these more adequate 
calculations; thus we hope it will survive for very strongly disordered pseudogaped su- 
perconductor as well. 

In Fig. [3n]we present the results of the computation of the T ~ Q differential Andreev 
conductivity based on Eas. (|200l20ip for several values of the contact transparency. Note 
that we give here the data for Andreev conductivity normalized to its value at zero 
voltage G(0), which scales itself as T^. Therefore the plots presented in FiglMl illustrate 
variation of the G{V) peak shape, but do not show the overall scale of G{V) as function 
of T. The gap function A(e) was computed using the T = Q limit of the solution to the 
modified MFA equation (|182|) for the "critical" M{uj), Eq.(|331). Voltage is normalized by 
the maximum value Aq of the gap function A(^). 

Note that the plots in Fig 135] illustrate the limit of a very large paring gap Ap 3> Aq. 
In the case of parity gap Ap comparable to the superconductive gap A(0), the point- 
contact conductance should contain two peaks: at eVi — A(0) due to Andreev processes, 
and at eV2 = Ap + A(0) due to single-electron tunnelling. These two voltage values 
behave differently with change of temperature and location of the point contact. The 
voltage Vi is determined by the collective superconductive gap Aq and, therefore, is 
position- independent but goes to zero as T — > Tc. In contrast, V2 fluctuates strongly 
between different contact locations due to the local paring ga p Ap but does not vanish 
as T approaches Tc from below. Recent experimental results 101| seem to support the 
above predictions. 



6.4- Spectral weight of high-frequency conductivity and superconducting density. 

The studics jl02| of high Tc superconductors has shown that the unconventional 
physics of these materials can be probed by the temperature dependence of the spec- 
tral weight K*°* (T) of the high-frequency conductivity defined by 

if *°*(T) = - / 3?a(c^, T)du; + p,(T) = K(T) + p,{T). (202) 
Jo 

Here a{uj,T) denotes the regular part of conductivity, whereas the contribution ps = 
—cS'il 5A.{uj = 0) of superconductive part (proportional to (5(w)) is explicitly separated 
in the last term of Eq. (|202p : the high-frequency cutoff ilmar is chosen a several times 
larger than Tc- 

In agreement with the data on conventional superconductors, BCS theory and its 
strong coupling generalizations predict that K*°*{T) does not change as temperature 
decreases below the superconducting transition temperature Tc- the appearance of Ps{T) 
is compensated by the decrease of the regular part. In contrast, for underdoped cuprates 
the regular part K{T) of the c-axis spectra l wei ght changes below Tc by one half of Ps{T) 
only, so that K*°\0) - K^°\Tc) = haJQ).^^ 

The observation jl02l | was explained by loffe and Millis )l03l | as being caused by 
formation of a pseudo-gap that survives up to a temperature TpQ well above Tc, and at 
the same time by the absence of the inter-layer coherence in the range Tc < T < TpQ due 
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to strong quantum and thermal fluctuations. The weight non-conservation occurs only 
for the c-axis cond uctivity, wh ereas the in-plane conductivity behaves in a usual way 
in all cuprates. |l02l | The theory |l03| attributes to the large phase fluctuations between 
weakly coupled layers and a tunneling nature of c-axis transport, in contrast with the 
smooth phase variations and the continuous electron motion in the planes. 

Wc expect that pseudogapped superconductivity formed near superconductor-insulator 
transition shows similar anomaly in the behavior of the spectral weight Ktot{T). In par- 
ticular, we show that inside the strong pseudo-gap region Tc ^ Ap the effect of the 
K*°*{T) non-conservation is even stronger: its variation AK*°*{T) is equal (within the 
mean held approximation) to ps{T), whereas the variations of the regular part of spectral 
weight AK{T) are smaller by ~ 1/Zcs factor. In this derivation we shall assume that 
pseudogap sets the largest scale in the problem, i.e. that the upper frequency cutoff 
^max satisfies Tc ^ ^^max ^ Ap. This assumption allows us to neglect completely sin- 
gle electron processes (which do not contribute for uj < Ap) and use pscudospin model 
for the computation of the conductivity in the entire frequency range. Furthermore, be- 
cause, as we shall show below the conductivity in the pscudospin model decreases fast at 
frequencies above Tc, we can replace the upper cutoff flmax oo when evaluating the 
spectral weight in this model. 

In order to compute the conductivity in the pscudospin model we need to include in 
the spin Hamiltonian (|145p the effect of the vector potential A: 

if = 5^ 20^1 - I ^ Af,, iS+S-e^'f"^ + h.c.) , (203) 

j ij 

where (pij = || J^l' Adr. Here i,j denote localized single-particle states (orbitals) with 
wavefunctions ipij{r). The Hamiltonian in the form (|203p is applicable provided that the 
typical distance Rq between centers of localization of the relevant cigenstatcs is much 
longer than the localization length Lioc (the same assumption is used everywhere in this 
section, see Eq. (|179p ). Then i,j can also be viewed as sites characterized by the position 
vectors ri,rj of centers of localization of the corresponding orbitals. The assumption 
Rq S> iioc implies that the theory is not sensitive to an exact definition of "center of 
localization" inside the localization radius. It is this assumption that allows to define the 
phase <f>ij in Ea. (|203p . Another assumption implied in (j203[) is that the magnetic field is 
weak enough so that its effect reduces to the introduction of the phase factors e**'^ and 
does not change the matrix elements Afy-. 

The electrical current across each link (ij) that corresponds to the Hamiltonian 
Ea. ipgg)) is 

4- = |a% *(5+5-e»^-^ - h.c). (204) 

The response of the superconductor to the slowly spatially ^arying electromagnetic field 
is fully characterized by Q{uj) = —cd\/6Ak=o{uj) = Qo + Q{uj) where we separated the 
constant part Qo and the frequency dependent one Q{io) that satisfies Q{ll: — > oo) ~ 0. 
At high frequencies the conductivity is purely imaginary with lma{Lu) = Qq/lo, expressing 
Im(T(w) through the real part by Kramers-Kronig relations one gets that K*°*{T) = Qq. 
The constant term Qo originates from the direct expansion of the expression for the 
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current (|204|) in (j)., 



K""{T) = Qo{T) = 




(205) 



where Qf/ ^ i^^^ ^■^■)- 

The equation (|205p imphcs that in the general case Qq cannot be temperature in- 
dependent because l^S^ SJ + /i.e.) is generally temperature dependent. Furthermore, 
one expects a strong temperature dependence to appear below Tc where order param- 
eter induces large spin-spin correlator. In the leading order in l/Z^ff one can replace 
{Sf Sj) « {S^) {Sj) and conclude that the spectral weight is constant above Tc but 
acquires temperature dependence below T^. 

K-totfrp^ 1 V- /^2e Y A,Aj-tanh/?£,tanh/3£j- 

K (T)-^1^9M^A-J-^^, — • (206) 

ij 

where = ^y £,f + Af. 

We now compute the leading correction in l/Zpff to (^S^Sj') above Tc , wc find 
that it does not result in the spectral weight temperature dependence in this tempera- 
ture range. This result does not preclude that temperature dependence appears in high 
orders in 1/Zeff but it implies that the temperature dependence is very small in this 
regime. Note that in this temperature range the superconducting response is absent, so 
temperature independence of the full weight implies temperature independence of its reg- 
ular part. To compute the leading contribution to (^S^S'^) it is sufficient to diagonalize 
two spin problem with the Hamiltonian 



H = ^,a^+(,a^~!^M,,{S+S-+S+S-), (207) 
which mixes states pairwise. We get 



{StS-+S+S-) = ^M., 



tanh(/36) - tanh(/?6) tanh(/36) + tanh(;36) 



4(6 - 6) 



4(6+6) 



Inserting this expression in (j205p and averaging over distribution of ^ , assuming that 
it is not correlated with the values of Mij and x^ we get 



Oo 



^Ah,] R 



R 



tanh(/36) - tanh(^6") 



2(6-6) 



The temperature enters only in the last factor in this product, differentiating it over the 
inverse temperature we find 



dR 
'dp 



^(6)^^(6) 
6-6 
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1 



1 



cosh^ /36 cosh^ ^6 



The integrand decreases exponentially fast for ^ ^ T , neglecting the contribution from 
these regions and assuming, as usual, that i>{^ < T) = v is constant at low energies we 
perform integration and get the result announced earlier: dR/dp ~ 0. 

We now discuss the regular and superconducting contributions to the spectral weight 
below Tc- The computation above demonstrates that leading terms in 1/Ze// involve 
only pairs of spins. This remains true below Tc as well and allows us to compute easily 
the regular part of the spectral weight in this temperature range. Wc begin by using the 
Kubo formula to express it through the spin-spin correlators 



dw 

— RcRij(uj) 



(208) 



where Ruiuj) is current-current correlator 



Rij(uj) = i 



dte" 



4(i)'4-(o) 



(209) 



evaluated in the absence of the external field: iij{t) = i{S^Sj — S^)t- The equation 
(|208p can be simplified by using the Kramers-Kronig relation for the correlator Rij{uj): 



^"(^) = ^ E (f lni%(0) 



(210) 



K{T < Tc) 



1 

8V 



E 



2e 



ij \ tanh(£j/r) - tanh(£j/r) 



ij \ tanh(£,/r) + tanh(£j/r) 



(211) 



gMi, 



One can immediately see that the regular part (|21ip contains an extra factor 
compared to the total spectral weight X*"* (T) in (|206p . Thus we conclude that variation 
of the regular spectral weight below Tc is smaller than total spectral weight: 



K{Tc) - K{T) gAU 



Tc 
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The last estimate in Ea. (|212p may be considered as a practical definition of Zcs- Indeed, 
in lattice spin models where each spin is coupled to Z^s other spins by a coupling constant 
gM the transition temperature can be estimated from gMZcS ^ Tc- As we argued 
above (see Eqs. p79ll8ip and the corresponding discussion in section WT\ . the effective 
coordination number Z^s is not too small (and superconductivity survives) in a wide 
region of localized single-particle states where 5l/Tc is larger than 1. Ea. (|212p implies 
that at Zea ^ 1 the temperature dependence of the total spectral weight at T < Tc is 
almost entirely related with the variation of the superconducting density: 



K^°\T) - K^°\Tc) 
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Figure 37: (Color online) Temperature dependence of regular spectral weight K{T) and full spectral 
weight K^°^{T) computed numerically (with 3D Gaussian AM of lattice size L = 20) in the region of 
pseudogap superconductivity at Ep = 7.9. We used Eas.l l211l l and II2U6II . correspondingly, and employed 
approximate expression for the gap function A(5, T) in the factorizcd form A(g, T) = Ac(§)\/1 — T/Tc, 
where Ac (5) is the critical-point solution shown in Fig. 1191 Black stars stand for the simple averaging, 
whereas blue circles correspond to typical averages. 



which is thus described by Ea. (j205[) . 

Wc now compute Ps{T) and K'^°*{T) as given by (|206p in the mean field approximation 
developed in this paper. We use modified MFA equation for the gap function in the form 
(|182p and replace xf^ under the sum by its average value ^R^ determined in Eq. (|179p in 
(pnS)) . The resuh is 



. HiSSjM r ..„h (213) 



which can be further simplified at T — > 0: 



2e2i?2 r Ke)rfe^2(^) ^ 2i^oe^^2 (214) 



h? e ' ' h? 

This expression allows one to extract the effective interaction range Rq from the data 
on the superconducting density in the combination with the gap A and DoS ^{e) found 
from the Andreev spectroscopy measurements (see section [6.3^ . 

In the above discussion we assumed that effective number of neighbors Z^s for the 
pseudospin model is large. It is not always the case because, as discussed in the end 
of section 16.11 the superconducting state survives in the strongly spatially fluctuating 
form in the regime Z^s < 1. Qualitative picture of the spectral weight temperature 
dependence in the pseudogap range of the Anderson model formulated in section [^?T] can 
be obtained from the numerical evaluation of the regular part K{T) given by (|21ip and 
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the total weight K*°\T) given by ((206)) . The results are shown in Fig 133 for Ep = 7.9. 
We emphasize that the expressions (j206l21ip give the leading orders in l/Zcn expansion 
only, which is the zero order for Eq. (|206p . and the first order for Eq. (|21ip . As one sees in 
Fig|571 both K{T) and K*°*'{T) change by a similar amounts, so in this example Zcs ^ 1. 
We expect that Fig[27] provides a qualitatively good description of the spectral weight 
behavior in pscudogapped superconductors. 



7. Summary of results and unsolved problems. 

We presented a generalization of the BCS-like theory of superconductivity that is 
appropriate for poor three-dimensional conductors in which Fermi level is located in the 
vicinity of the mobility edge. The developed theory neglects the effects of Coulomb re- 
pulsion. The reasons why this approximation is appropriate to describe many disordered 
films were discussed in section 11.21 

The model with BCS-likc attraction and no Coulomb repulsion solved in the bu lk of 
the paper can be also real ized by ultra-cold fermionic atomic gases in optical traps [6l. ll04j 
with controlled disorder jl05| . The enhancement of by disorder predicted theoreti- 
cally might turn out to be a useful tool to observe BCS state in these systems because 
without such enhancement Tc is often too low to be reached experimentally for the small 
values of the interaction constant A. Such experiments may be also important for the 
understanding of superconductor-insulator transition in general. 

The control parameter of the model is position of the Fermi level Ep with respect to 
the mobility edge Ec- The major new (with respect to the old works 0, d, 3) ingredient 
of the developed theory is the full account of the critical and weakly off-critical wave- 
function's fractal! ty. We identified three qualitatively different regimes: the hard-gap 
insulator realized at large disorder, the pseudogap superconductor that appears when 
the phase correlations develop between well localized Cooper pairs and the "fractal" 
superconductor that appears against the background of fractal single particle states. 
Upon a further decrease of disorder the fractal superconductor smoothly crosses over to 
a conventional one. 

For the insulating state, our main result is the computation of the activation en- 
ergy Tj in the Arrhenius temperature dependence of resistivity, lni?(r) ~ Tj/T, at low 
temperatures. The hard-gap character of this insulating state is due to the local attrac- 
tive electron-electron interaction which leads to formation of localized Cooper pairs with 
nonzero binding energies Ap-*, specific for each j-th localized orbital. The probability 
distribution of these energies is similar to the distribution of inverse participation ratios 
for localized eigenstates of the 3D Anderson model; it drops exponentially fast at low 
energies, imitating a hard gap. The estimate for the activation energy T/ as a function 
of Ep — Ec is in a reasonable agreement with experimental data (iij . 

The fractal superconducting state (realized for Ep very close to Ec) is characterized 
by the following features: 

1. In the weak coupling limit, A ^ 1, the transition temperature, Tc, becomes 
a power law function of the interaction constant, A. This leads to the parametric 
enhancement of Tc with respect to its value deep in the metal state. This conclusion 
is different from previous works which assumed the validity of the "Anderson 

theorem" in the Ep region very close to mobility edge. The power law exponent in 
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the Tc{X) oc A^/*^^ ''^^ dependence is determined by the fractal dimension d2 ~ 1-3 
of the critical eigenfunctions. 

2. Strong local fluctuations of the pairing amplitude A(r) coexist with a unique and 
well-defined critical temperature Tc below which a macroscopically coherent state 
appears. 

3. The thermal Ginzburg parameter Gi of this transition is a universal quantity of 
the order of unity which is independent of the interaction constant A, the same 
holds for the mesoscopic Ginzburg number Gi^. This implies that thermal and 
mesoscopic fluctuations may change the pre-factor in the power law Tc{X) by a 
factor ~ 0(1) but cannot change the power-law itself. 

4. Local single-particle density of states z^(e, r) fluctuates strongly in real space; these 
fluctuations lead to a strong random asymmetry of the tunneling conductance, 
GT{V,r) ^ GTi-V,r). 

5. The value of the supcrfluid response, Ea. (jl76p . coincides with result of (Equa- 
tion (4.1) in this paper), provided that one inserts in it the correct values of Tc and 
Aq given above. 

The pscudogap superconductor is predicted to occur when the Fermi energy Ep is 
deep inside the localized band of single-particle states. Its major feature is the presence 
of two distinct energy scales, both originating from the Cooper attraction between weakly 
localized electrons: the collective superconductive gap A(0) and the local binding energy 
of a Cooper pair Ap. The most unusual behavior is expected to occur in the regime of a 
strong pseudo-gap Ap ^ (A(0), Tc). The very existence of superconducting correlations 
and a nonzero Tc in this regime is unexpected because it is characterized by a typical 
level spacing 5l = [vqL^ioc)'^ which is larger than Tc and A(0): Sp/Tc (x(Ap/Tc)^^'^^ 
(|81l82p . The appearance of this regime was not expected in previous studies which 
concluded that the superconducting state is stable only up to dp l£ Tc- The perseverance 
of superconducting coherence much deeper in the localized region than was expected pre- 
viously is the result of the enhancement of the correlations between the wave-functions 
intensities [t^I, which occurs due to the Mott's mechanism of resonant mixing of localized 
states. The key features of the pseudo-gap superconductor (in addition to the features 
of fractal superconductor listed above) are the following: 

1. Insulating behavior of the resistivity in a wide range of temperatures above Tc. 

2. Formation of almost hard gap without coherence peaks in the density of states 
above Tc with the gap value of the gap that fluctuates signiflcantly from point to 
point. 

3. Growth of coherence peaks below a global Tc. The magnitude of the coherence 
peaks (proportional to the local value of the superconducting order parameter) 
fluctuates from point to point, these fluctuations become very large close to the 
superconductor- insulator transition. 

4. Two-peak feature in differential conductance at moderate transmission probabili- 
ties as measured by the Andreev point-contact spectroscopy below Tc- The lower 
peak voltage Vi{T,r) is expected to be r -independent, but vanishing as temper- 
ature approaches Tc from below. In contrast, the higher peak voltage V2{T,r) is 
expected to be almost T-independent, but strongly r-dependent. 
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5. Strong temperature dependence of the total spectral weight K*°*{T) of the high- 
frequency conductivity which shows with a loss of major part of K^°^{T) at tem- 
peratures T < Ap, and its re-appcarancc below Tc. The regular part K{T) of the 
spectral weight is smaller than the superconducting response ps, sec Ea. (j212p . 

We now briefly discuss possible extensions of the developed theory and open questions. 
We begin with purely theoretical questions. 

An obvious extension of the present theory would be a consistent account of Coulomb 
interaction effects in fractal superconductors near the 3D mobility edge. Here one should 
distinguish the effects of short range part of the Coulomb repulsion that competes with 
the phonon attraction and the long range part which might become very important in 
the insulator resulting in the formation of the Cooper pair glass and qualitatively new 
physics. Decreasing the dielectric constant in a real material would change it from the 
fractal superconductor discussed in this paper to the material in which superconductivity 
is suppressed by Coulomb interaction. Furthermore, increasing only the short range part 
of the Coulomb interaction leads to the 'fermionic' mechanism while increase in the long 
range part leads to the Coulomb blockade with completely different properties. The 
transition between all these regimes are currently not understood, neither theoretically 
nor experimentally. It seems likely that at least one of these crossovers was observed as 
the 'region of poor scaling' in [62| . 

The character of the quantum critical point that separates the fractal superconductiv- 
ity and the insulator at T = deserves a further study. As discussed in section IH] on both 
sides of this transition electrons are bound in the localized pairs, and the transition itself 
consists in development of phase correlations, well described in terms of the XY ordering 
of the Anderson pseudo-spins Sj. This scenario is different from both "bosonic" and the 
"fermionic" mechanisms of superconductor-insulator transition, so we suggest the name 
"pseudo-spin" mechanism. The cavity approach to the study of such transitions was 
developed recently in [o^ ; the main qualitative conclusion of this study is self-organized 
inhomogeneity of the resulting superconductor in which phase correlations are dominated 
by the rare, almost one-dimensional paths. This conclusion is supported by the results 
of the virial expansion method shown in Fig. [29] that point out to a percolation-like 
transition between the pseudo-gap superconductor and an insulator. Another feature of 
this solution is a very rapid decrease in the transition temperature beyond certain value 
of disorder which might be seen experimentally as the apparent existence of the lowest 
nonzero Tc ~ 0.5K found in amorphous InO^^ system [42l |. 

We now turn to experimental findings that lack (partially or entirely) theoretical 
explanation. 

As discussed in section 16.2.31 the quantitative description of point contact tunneling 
requires the theory of collective modes in the fractal superconducting state. In the 
absence of such theory the spectacular, nearly rectangular shape of the local tunneling 
conductance G{V) at low temperatures [i^l cannot be explained quantitatively. 

The theory of the hard-gap insulating state presented in section |3] neglects the trans- 
port by Cooper pairs and takes into account only single electron transport. It is therefore 
limited to the region relatively far from the transition such as studied in Ref. [i^ . Theory 
of the incoherent transport by Cooper pairs close to the transition was discussed in the 
recent papers [H, Il06l | , the main conclusion of these studies is that very close to the 
transition the behavior of Cooper pairs is controlled by the many body mobility edge 
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which becomes zero exactly at the transition. Away from the transition the theory pre- 
dicts activated behavior at the lowest temperatures with the gap that increases fast away 
with disorder and becomes infinite a short distance away from the transition signalling 
the absence of incoherent transport by C ooper pairs [93 1 . 

Recently the experimental paper |107l | reported the anomalous size-dependence of the 
superconductor- insulator transition for extremely wide range of sizes between 1 and 150 
microns. In particular, the workers observed a slow but significant dependence of the 
activation energy Tq{L) on the system size L (the distance between the metal contacts) 
on the insulating side of the transition. In one case Tq{L) decreased from 13.5 K to 9.6 
K while L was changed from 145 to 12 microns at a fixed width. Although the theory of 
the the incoherent transport by Cooper pairs close to the transition 93| predicts a strong 
size dependence at mesoscopic scales due to inhomogeneities, the size dependence at the 
huge scales observed experimentally arc very difficult to explain. The explanation of this 
effect seems currently to be beyond the reach of a theory and presents a real challenge. 

A different set of challenges is presented by the superconductor-insulator transition in- 
duced by magnetic field for strongly disordered (almost insulating) superconducting sam- 
ples. As one would expect, the critical field. He (which corresponds to superconductor- 
insulator transition at T = 0) is a strong function of disorder, for instance, it varies 
over two orders of magnitude in the experiments [Hsj . Surprisingly, there are indications 
that phenomen ology o f the transitions driven by disorder and by magnetic field look 
differently [46, 6l[ 44, 45, 62, 65 1. Namely, the critical point in disorder is character- 
ized by the activated behavior of the resistance, R ~ exp(To/T), with a large gap Tq 
even very close to the superconductor-insulator transition. In contrast, activated behav- 
ior, R ^ cxp{Tj{H)/T), observed at fields above critical He is characterized by a small 
Ti{H) which extrapolates to zero at He ■ In other words, the region in {H,T) plane 
where transport is characterized by activated behavior with a small gap has a peculiar 
wedge like shape as sketched in Figured]). Generally, one expects that very close to the 
critical line the transport is dominated by Cooper pairs; this expectation was experimen- 
tally confirmed for ultrathin Bi films [68[. In this case the experimental phase diagram 
implies that the regime of Cooper pair dominated transport is narrow in disorder but 
wide in field. Qualitatively this is likely to be due to the fact that the effect of a small 
magnetic field is limited to the generation of local phase differences which have a large 
effect of the superconducting state but does not affect the incoherent transport of the 
pairs. However, there is currently no consistent theory of this effect. On the experimental 
side, it would be important to verify that the transport in this regime is due to Cooper 
pairs for the InO and TiN films. 

Furthermore, strongly disordered superconductors show puzzling behavior in very 
large magnetic fields: their resistivity become temperature independent and approaches 
the quantum limit, h/e'^. The field scale, i?p, at which this behavior sets in is not 
sensitive to the disorder in contrast to the critical field. He- The second field scale appears 
first in the fractal superconductor: it is the field that destroys local Cooper pairs and 
suppresses the hard gap Ap. In this regime, the parametric difference between Hp and 
He is due to the low fractal dimension o?2 ~ 1-3 of a single-particle wave- function that 
suppresses the orbital effect of magnetic field on the pairing of two electrons localized on 
this eigenstate, similarly to a very thin wire. In contrast, long-range coherence requires 
the existence of large loops involving many localized Cooper pairs. The typical size of 
these loops is larger than the localization length Lioe, leading to relatively small magnetic 
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field scale He < ^o/^Lc ^ ^p- ^^e pseudogap superconductor the second field scale is 
associated with the suppression of a large gap, Ap, responsible for Cooper pair formation 
whereas the critical field is associated with the frustration of a weak pseudospin coupling. 
Thus, one expects even larger difference between the field scales in this regime. However, 
the universal value of the conductance in the regime B ^ Hp is not expected theoretically 
and remains mysterious. 

Ncrnst effect can be potentially a very sensitive probe of the nature of the super- 
conducting state. |l08l | These experiments were performed recently on InO^: samples and 
show |109| Nernst signal which scales as N{T) oc T"" , with the "Nernst exponent" 
n « 7.5 at low magnetic fields. Application of the conventional th eory 18911 of the super- 
conducting fluctuations to the Nernst effect in 2D superconductors llOl . lllj gives Nernst 
signal that scales as l/TlnT at high temp eratures. This behavior was indeed observed 
in conventional superconductor NbSi [112|. Because Nernst effect requires a motion of 
the electrical charge around the plaquette, the temperature dependence is expected to 
be different in discrete systems where motion around the minimal plaquette involve a 
large number of hop s and each hop implies an extra power of 1 /T in the high tempera- 
ture expansion [ll3|. For instance, because minimal plaquette on the hexagonal lattice 
contains six sites, one expects that N{T) oc in this case. The striking difference 
between NbSi and InO^; behavior observed in [lOOj on very low Tc (very disordered) sam- 
ples indicates, in our opinion, the importance of Cooper pair hoping between localized 
sites. This behavior should get less pronounced further away from the transition, thus 
we expect that slightly less disordered samples will show more 'conventional' exponent 
in Nernst effect. 

A number of papers noted the apparent similarity between phenomenology of disor- 
dered film s of c onventional superconductors discussed in this paper and that of high Tc 
oxides. 53|, |114| Very briefly, the transport and magnetic measurements show the forma- 
tion of ps eudo gap at temperatures much higher than superconducting Tc in undcrdoped 
materials. |ll5l | while STM measurements display highly inhomogencous order parame- 
ter ( as ineasured by coh erence peaks) combined with modestly homogeneous tunneling 



gap. |116l . IllTl . IllSl Ill9l | Furthemore, a number of indirect e vidences points out to the 



;es p 

superconducting nature of the pseudogap in these materials [120|, 12ll . Il22| , similar to 
the situation discussed in this paper. The crucial difference between s-wave and d-wave 
pairing characterizing high Tc oxides is that, in contrast to s-wave pairing, modest elas- 
tic scattering with mean free path of the order of superconducting coherence length 
suppresses the d-wavc superconductivity. Thus, localization of the wave functions is in- 
compatible with generic d-wavc superconductivity. This conclusion might need to be 
revised for the special case of superconductivity which is due to pairing of electrons in a 
small area of the full Brillouin zone which might be the case of cuprates 123l 124l |. In 
this case, the electrons that are responsible for the pairing belong to two well separated 
patches on the Fermi surface with each patch characterized by a small momentum po- 
The density of states in these patches is large which makes scattering of these electrons 
strong. Furthermore, because most impurities are located far from the conducting copper 
oxide planes, it is likely that the elastic scattering, though strong, does not mix differ- 
ent patches and thus does not inhibit d-wave pairing. In this situation, single electron 
states inside each patch may become very similar to the localized states discussed in 
this work making the developed theory qualitatively correct. This would explain the 
main phenomenological features mentioned above: formation of pseudogap far above T^ 
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highly inhomogeneous order parameter observed in tunnehng data, insulating behavior 
of LaSrCuO superconductors in high magnetic fields, etc. 

Finally, it is not clear why the phenomenology displayed by Josephson junction 
networks (oi [soj was not observed in any disordered film. In particular. Josephson junction 
networks close to superconductor-insulator transition are characterized by the appear- 
ance of a wide region of magnetic fields where resistivity has no temperature dependence 
but varies by many orders of magnitude as a function of magnetic field. This is in a sharp 
contrast with the films that show more or less good crossing point of R{B) isotherms. 
There is very little doubt that Josephson networks are described by model Hamiltonian 
so one has to conclude that it is not the appropriate model for most films. One 
possible source of difference discussed above might be that many films are characterized 
by a large value of the dielectric constant, while another reason, relevant for the films 
with small dielectric constant, might be the capacitance matrix is not dominated by 
nearest neighbors in these films. Nearest neighbor capacitance matrix translates into the 
Coulomb interaction that decreases only logarithmically with distance in the Josephson 
network as opposed to 1/r expected in a film. The physical result might be the forma- 
tion of the Cooper pair glass in the Josephson network but not in the disordered film or 
vice versa. If true, one expects that Josephson networks with the Coulomb interaction 
screened by the ground plate might display properties more similar to those of the films. 
Another difference might be due to the fiuxes produced by magnetic field are completely 
random in a disordered films but are (almost) identical on the plaquettes of the Josephson 
lattice. In this case one expects to see film behavior in randomized Josephson network. 
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Appendix A. Virial expansion in pseudospin subspace 

The two-spin susceptibility — + Xi2 ^ is computed as explained in section 
14.41 and can be easily written analytically even for nonzero Hartree terms. It is a sum of 
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two terms, one of which does not exphcitly depend on Jj^: 
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Here we present the details of the formalism for computing the three-spin suscepti- 
bility. It is developed with the goal to use for numerical computation rather than being 
a full-analytical approach. The even Hilbert subspace of the problem of three sites with 
2 or electron on each site (equivalent to the spin-i problem with three coupled spins) 
consists of 8 basic states classified by the total spin. The corresponding Hamiltonian has 
a block structure corresponding to the total spin |, 0, ^ 
just numbers: 

±(Ci + e2 + e3)-^(ji'2 

The last two blocks are 3x3 matrices: 



i. The first two blocks are 



Jll, + J, 
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where = S.b + — ■ Such structure suggests that the eigenvalues of the 3-spin 
problem are also grouped: 



x(o) 



±(e. + + 60 - ^(^2 + jL + jL) 



(A.6) 

are found from the solution of the two 



a1°1 5 /) = 0, and dct(B - a[j°^ g /) = 0. 



and two groups of eigenvalues X^'^ 

cubic characteristic equations det(A — /vg 4 , 

One can find the solutions to these equations either numerically or from the Cartan 
formula for the roots of cubic equations, the former being more efficient in the numerical 
calculations. 

The 8-eigenvectors Xi = (1,0, 0, 0, 0, 0, 0, 0) and X2 = (0, 0, 0, 0, 0, 0, 0, 1) are simple: 
they contain unity in the first and eight row respectively, and all the other elements being 
zero; the structure of the nontrivial eigenvectors of the 3 x Smatrices is as follows: 



^1=3,4,5 = (0 , Uq , , Ifa , , , , 0) 
^1=6,7,8 = (0,0, 0,0<,<,<,0) 



(A.7) 



where {ua,Va,Wa) is the normalized 3-eigenvector of the matrix A that corresponds to 
the eigenvalue A^^la 4 5 and u'^,v'^,w'^ is the normalized 3-eigenvector of the matrix B 

that corresponds to the eigenvalue A^^lg 7 g . 

The perturbation matrix in the same basis takes the form: 
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With the eigenvalues Xa ^ and the eigenvectors Xa (a = 1. 
-fa from Eq. ([15|) (the field A cancels out): 

{up + up + wpY 
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8) obtained, one can compute 
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and then compute xijk from Ea. (jl54p . Finally, subtracting the proper combinations 
of the two- and one- spin susceptibilities one obtains the virial coefficient X3(T) from 
Eq.dinil). 
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Appendix B. Virial expansion including single-occupied states 



The Hilbert space we considered so far corresponded to each orbital i,j,k... occupied 
either by two electrons or empty. Now let us consider the same Hamiltonian, Ea. (jl45|) . 
but extend the Hilbert space to include possible single-electron occupancies of all orbitals. 
To avoid unnecessary complications we omit here the Hartree terms setting g|| = in 
Eq.dHSl). 

Appendix B.l. One- orbital problem 

The full Hilbert space in this case consists of 

|0>, |2), It), U). 

The first two states contain even number of electrons (the "even" states) and the last two 
(the "odd" states) contain one electron each. Each "odd" state has extra energy Ap. 

It is important that the spin operator acting on the "odd" states | f) or | l) gives 
zero: 

5+/-|;) =5+/-|t) =0. (B.l) 
The energies of the four states are (we count all energies with respect to — e = — ^ + Ap): 

-e, 2C-2Ap = ^-£ = Ap, Ap. 

Ea. (jB.l|) suggests that the "partial susceptibilities" 7odd given by Ea. (|155|) are zero for 
all the odd states. This is because the Hamiltonian Ea. (|145|) docs not couple the even 
and odd sectors of the Hilbert space. 

Thus all one has to do to compute the one-site susceptibility is to add 2e~^^l^ to 
the partition function Zq. The result is: 

,... = (^)^„,A,,, ,B., 



Note that there is a nice property: 

F(£,oo)=tanh^-j , ^^(e, 0) = tanh ^— j . (B.3) 

Appendix B.2. Two-orbital problem 

The full Hilbert space for the 2-orbital problem consists of the even sector: 

^i = |0,0), V^2 = |2,2), V3 = |2,0), ^4 = |0,2) 

an the odd sector: 

V's.e = |a,0), i/jT^s = \0,a), ^^940 = I", 2), 

■011,12 = |2,a), V13-16 = !«,«')> 

where a or t- 
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Let us first find the eigenvalues in the odd sector (counted from — (ei + £2)): 



Xfl = 6-(£i+e2) = -ei + A^'\ 
>^9A0 = ei+2(6-A^'')-(ei+£2) = e2 + A«, 
^n]i2 - 2(ei-Ai,'))+6-(ei+£2) = ei + A^''. 
In all the above formulae 

e.=6-A«, (B.4) 

where is a single-particle energy (eigenvalues of the state '0i(r)) measured from the 
Fermi-energy, and A^p is a pseudo-gap: 

r 

Note that Ea. ()B.ip guarantees that the eigenvalues of the odd sector do not have any 
My-dependent renormalization. However, not all the partial susceptibilities 7odd are 
equal to zero. The reason is that the source terms A S'^ + A* S~ contain only one 
5"*"'"^ operator whereas the number of sites is two. This operator may act on the site 
in an even state and then the result of action of the source term on a semi-odd state 
|even, odd) will be non-zero. Thus we expect that the structure constants for the states 
ipi, {i = 5, 6, ...12) are non-zero. A simple calculation shows that 

1 

75 = 76 = -79 = -710 = -7^, 

2e2 

1 

77 = 78 = -711 = -712 = ~-z — ■ 

2ei 

Now we are in a position to give an exact result for the 2-orbital susceptibility: 



(2) _ 2 _E+_E_ ' 2_E+ ^ 2E_ ' 2ei ^ 2s2 

~ cch(^,^) +cosh(f ) e-^pVT + cosh(f) e-<VT + e-A<,^)/Tg-Aj?>/T 

(B.6) 



where we denoted 



ssh(a;, y) ~ sinh (x) sinh (y) , sch(x, y) ~ sinh (a;) cosh (y) , cch(a;, y) = cosh (x) cosh (y) 

If one neglects the renormalization of eigenvalues by interaction and sets i?+ = ei, 
= 62, the Ea. (jB.6p is reduced to a simpler form 



X?^ = -^n=„AW)F(e2,A(?)) + 

Z£ie2 

F(ei,AW) , F(e2,A(?)) 



2£i 2£2 
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where F{e,Ap) is defined in Ea. ()B.2p . In this approximation the fuh susceptibihty has 
the same form as Ea. (jl58p after one replaces tanh(£/T) by the function F{e,Ap). 

However, when the renormahzation of energy levels is taken into account, such simple 
replacement tanh(i?/r)— > F{E,Ap) does not work any more. 

Appendix B.3. Three-orbital problem 

In an absolutely similar way one can compute the exact susceptibility of a three-site 
problem, considering all 64 states of which 8 are even and 56 odd states are grouped as 
follows: 
group I: 

I two — site even state) (K) |a) 

with a =t (i) corresponding to the first, second or the third site. There are 4 x 2 x 3 = 24 
such states. 
group II: 

|one — site even state) ® \a) (g) \a'}, 

with an even state on the first, second, or the third site. There are 2 x 4 x 3 = 24 such 

states. 

group III: 

\a, a' a"). 

There are 8 such states. 

The eigenvalues Xa^ can be found exactly like in the previous subsection; they are 
equal to (relative to the vacuum state — (ei-l-£2+£3)) the eigenvalues of the corresponding 
even state plus the sum of Ap ' of the singly occupied sites. For instance: 

|0, 2, a) ^ A(") = -sgn(£i - £2) - £2)^ + | J12P + A^^), 

|2,a,a')->AW=ei+A(2)+A(3), 

|a,a',a")->A(0)=AW+A(2)+A(3). 

The three-orbital susceptibility xfl 3 be calculated using the general Ea. ()154p . in 
which the partition function Zq = Zq(1, 2, 3) is equal to 

Zo(l, 2, 3) = + 8 cosh{E^^'^^) cosh(^i'''^) e'^p' 

+8 cosh(^*!"''') cosh(^i^'^') e-^p^ + 8 cosh(^^^'^)) cosh(^i^'^^) e^'^p" (B.8) 
+ 8cosh(£i)e e 4- 8 cosh(£2) e ^p e ^p -F8cosh(£3)e ^p e ^p 
+ 8e p e p e p . (-d-9) 

In Ea. (|B.8p each line corresponds to the above group of states, E± is given by 
Eq.dinni), and E,A = E/T,A/T. The partition function of the even states Zg™"^ = 

Seven states e~^° '/"^ has been calculated above for the 3-spin case. 
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The numerator of Eq. (|154p 0(1,2,3) = — ^ ^°'^"^7q is found to be the following: 

G(l, 2, 3) = G<=™"(1, 2, 3) + 4 G12 e-^p' + 4 G13 er^p' + 4 G23 e"^^" 
£1 £2 £3 



where 



G„- = sinh(^l''^'h sinh(^'^-^^)+ (B.ll) 



sinh(i;l*'^'' ) , - sinhfS^''^' 1 ,-u,\^ 
^,.\ ' cosh S^''^' + cosh S^'^' 



is (up to a constant factor) the numerator in Eq. (|154p for the 2-spin problem. Note also 
that sinh(£i)/ei is (again up to a constant factor) the numerator of the one-spin problem. 
The quantity G°™"(1, 2,3) = — X^cvon states e~'^° '/^ 7^ has been calculated above for the 
3-spin problem. 

Using Eas. (|B.8p . (|B.10p one finds the 3-orbital susceptibility exactly: 

(3) G{1,2,3) 

- Zo(l,2,3)- ^^■^^> 
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